Page MenuHomeHEPForge

No OneTemporary

diff --git a/Models/General/BSMModel.cc b/Models/General/BSMModel.cc
--- a/Models/General/BSMModel.cc
+++ b/Models/General/BSMModel.cc
@@ -1,365 +1,365 @@
// -*- C++ -*-
//
// This is the implementation of the non-inlined, non-templated member
// functions of the BSMModel class.
//
#include "BSMModel.h"
#include "ThePEG/Interface/ClassDocumentation.h"
#include "ThePEG/Interface/Parameter.h"
#include "ThePEG/Interface/Switch.h"
#include "ThePEG/EventRecord/Particle.h"
#include "ThePEG/Repository/UseRandom.h"
#include "ThePEG/Repository/EventGenerator.h"
#include "ThePEG/Utilities/DescribeClass.h"
#include "ThePEG/Repository/Repository.h"
#include "ThePEG/Utilities/StringUtils.h"
#include "ThePEG/Persistency/PersistentOStream.h"
#include "ThePEG/Persistency/PersistentIStream.h"
#include "ThePEG/PDT/MassGenerator.h"
#include "ThePEG/PDT/WidthGenerator.h"
#include "ThePEG/PDT/DecayMode.h"
using namespace Herwig;
BSMModel::BSMModel() : decayFile_(), readDecays_(true),
topModesFromFile_(false),
tolerance_(1e-6)
{}
void BSMModel::persistentOutput(PersistentOStream & os) const {
os << decayFile_ << topModesFromFile_ << tolerance_;
}
void BSMModel::persistentInput(PersistentIStream & is, int) {
is >> decayFile_ >> topModesFromFile_ >> tolerance_;
}
DescribeAbstractClass<BSMModel,Herwig::StandardModel>
- describeHerwigBSMModel("Herwig::BSMModel", "");
+describeHerwigBSMModel("Herwig::BSMModel", "Herwig.so");
void BSMModel::Init() {
static ClassDocumentation<BSMModel> documentation
("The BSMModel class provides a base class for BSM models including the"
" features to read decays in the SLHA format");
static Parameter<BSMModel,string> interfaceDecayFileName
("DecayFileName",
"Name of the file from which to read decays in the SLHA format",
&BSMModel::decayFile_, "",
false, false);
static Switch<BSMModel,bool> interfaceTopModes
("TopModes",
"Whether ro use the Herwig++ SM top decays or those from the SLHA file",
&BSMModel::topModesFromFile_, false, false, false);
static SwitchOption interfaceTopModesFile
(interfaceTopModes,
"File",
"Take the modes from the files",
true);
static SwitchOption interfaceTopModesHerwig
(interfaceTopModes,
"Herwig",
"Use the SM ones", false);
static Parameter<BSMModel,double> interfaceBRTolerance
("BRTolerance",
"Tolerance for the sum of branching ratios to be difference from one.",
&BSMModel::tolerance_, 1e-6, 1e-8, 0.01,
false, false, Interface::limited);
}
void BSMModel::doinit() {
StandardModel::doinit();
// check if need to read decays
if(decayFile()==""||!readDecays_) return;
decayRead();
}
void BSMModel::decayRead() {
// read decays
CFileLineReader cfile;
cfile.open(decayFile_);
if( !cfile ) throw SetupException()
<< "BSMModel::doinit - An error occurred in opening the "
<< "decay file \"" << decayFile_ << "\"."
<< Exception::runerror;
//Before reading the spectrum/decay files the SM higgs
//decay modes, mass and width generators need to be turned off.
PDPtr h0 = getParticleData(ParticleID::h0);
h0->widthGenerator(WidthGeneratorPtr());
h0->massGenerator(MassGenPtr());
h0->width(ZERO);
h0->stable(true);
DecaySet::const_iterator dit = h0->decayModes().begin();
DecaySet::const_iterator dend = h0->decayModes().end();
for( ; dit != dend; ++dit ) {
generator()->preinitInterface(*dit, "BranchingRatio", "set", "0.");
generator()->preinitInterface(*dit, "OnOff", "set", "Off");
}
// if taking the top modes from the file
// delete the SM stuff
if(topModesFromFile_) {
PDPtr top = getParticleData(ParticleID::t);
top->widthGenerator(WidthGeneratorPtr());
top->massGenerator(MassGenPtr());
DecaySet::const_iterator dit = top->decayModes().begin();
DecaySet::const_iterator dend = top->decayModes().end();
for( ; dit != dend; ++dit ) {
generator()->preinitInterface(*dit, "BranchingRatio", "set", "0.");
generator()->preinitInterface(*dit, "OnOff", "set", "Off");
}
}
// read first line and check if this is a Les Houches event file
cfile.readline();
bool lesHouches = cfile.find("<LesHouchesEvents");
bool reading = !lesHouches;
if(lesHouches) cfile.readline();
// function pointer for putting all characters to lower case.
int (*pf)(int) = tolower;
while (true) {
string line = cfile.getline();
// check for start of slha block in SLHA files
if(lesHouches && !reading) {
if(line.find("<slha")==0) reading = true;
if(!cfile.readline()) break;
continue;
}
// ignore comment lines
if(line[0] == '#') {
if(!cfile.readline()) break;
continue;
}
// make everything lower case
transform(line.begin(), line.end(), line.begin(), pf);
// start of a block
if(line.find("decay") == 0) {
readDecay(cfile, line);
continue;
}
else if( lesHouches && line.find("</slha") == 0 ) {
reading = false;
break;
}
if(!cfile.readline()) break;
}
}
void BSMModel::readDecay(CFileLineReader & cfile,
string decay) const{
// extract parent PDG code and width
long parent(0);
Energy width(ZERO);
istringstream iss(decay);
string dummy;
iss >> dummy >> parent >> iunit(width, GeV);
- PDPtr inpart = getParticleData(parent);
+ PDPtr inpart = getBSMParticleData(parent);
if(!topModesFromFile_&&abs(parent)==ParticleID::t) {
cfile.readline();
return;
}
if(!inpart) throw SetupException()
<< "BSMModel::readDecay() - "
<< "A ParticleData object with the PDG code "
<< parent << " does not exist. "
<< Exception::runerror;
inpart->width(width);
if( width > ZERO ) inpart->cTau(hbarc/width);
inpart->widthCut(5.*width);
Energy inMass = inpart->mass();
string prefix(inpart->name() + "->");
double brsum(0.);
unsigned int nmode = 0;
while(cfile.readline()) {
string line = cfile.getline();
// skip comments
if(line[0] == '#') continue;
// reached the end
if( line[0] == 'B' || line[0] == 'b' ||
line[0] == 'D' || line[0] == 'd' ||
line[0] == '<' ) {
cfile.resetline();
break;
}
// read the mode
// get the branching ratio and no of decay products
istringstream is(line);
double brat(0.);
unsigned int nda(0),npr(0);
is >> brat >> nda;
vector<tcPDPtr> products,bosons;
Energy mout(ZERO),moutnoWZ(ZERO);
string tag = prefix;
int charge = -inpart->iCharge();
while( true ) {
long t;
is >> t;
if( is.fail() ) break;
if( t == abs(parent) ) {
throw SetupException()
<< "An error occurred while read a decay of the "
<< inpart->PDGName() << ". One of its products has the same PDG code "
<< "as the parent particle. Please check the SLHA file.\n"
<< Exception::runerror;
}
- tcPDPtr p = getParticleData(t);
+ tcPDPtr p = getBSMParticleData(t);
charge += p->iCharge();
if( !p ) {
throw SetupException()
<< "BSMModel::readDecay() - An unknown PDG code has been encounterd "
<< "while reading a decay mode. ID: " << t
<< Exception::runerror;
}
++npr;
tag += p->name() + ",";
Energy mass = p->mass();
mout += mass;
if(abs(p->id())==ParticleID::Wplus||p->id()==ParticleID::Z0) {
bosons.push_back(p);
}
else {
products.push_back(p);
moutnoWZ += mass;
}
}
if( npr != nda ) {
throw SetupException()
<< "BSMModel::readDecay - While reading a decay of the "
<< inpart->PDGName() << " from an SLHA file, an inconsistency "
<< "between the number of decay products and the value in "
<< "the 'NDA' column was found. Please check if the spectrum "
<< "file is correct.\n"
<< Exception::warning;
}
if( npr > 1 ) {
tag.replace(tag.size() - 1, 1, ";");
if(charge!=0) {
cerr << "BSMModel::readDecay() "
<< "Decay mode " << tag << " read from SLHA file does not conserve charge,"
<< "\nare you really sure you want to do this?\n";
}
++nmode;
// normal option
if(mout<=inMass) {
inpart->stable(false);
brsum += brat;
createDecayMode(tag, brat);
}
// no possible off-shell gauge bosons throw it away
else if(bosons.empty() || bosons.size()>2 ||
moutnoWZ>=inMass) {
cerr << "BSMModel::readDecay() "
<< "The decay " << tag << " cannot proceed for on-shell "
<< "particles, skipping it.\n";
}
else {
Energy maxMass = inMass - moutnoWZ;
string newTag = prefix;
for(unsigned int ix=0;ix<products.size();++ix)
newTag += products[ix]->name() + ",";
if(bosons.size()==1) {
cerr << "BSMModel::readDecay() "
<< "The decay " << tag << " cannot proceed for on-shell\n"
<< "particles, replacing gauge boson with its decay products\n";
vector<pair<double,string> > modes =
createWZDecayModes(newTag,brat,bosons[0],maxMass);
for(unsigned int ix=0;ix<modes.size();++ix) {
modes[ix].second.replace(modes[ix].second.size() - 1, 1, ";");
createDecayMode(modes[ix].second,modes[ix].first);
brsum += modes[ix].first;
}
}
else if(bosons.size()==2) {
bool identical = bosons[0]->id()==bosons[1]->id();
if(maxMass>bosons[0]->mass()&&maxMass>bosons[1]->mass()) {
cerr << "BSMModel::readDecay() "
<< "The decay " << tag << " cannot proceed for on-shell\n"
<< "particles, replacing one of the gauge bosons"
<< " with its decay products\n";
unsigned int imax = identical ? 1 : 2;
if(imax==2) brat *= 0.5;
for(unsigned int ix=0;ix<imax;++ix) {
string newTag2 = newTag+bosons[ix]->name()+',';
unsigned int iother = ix==0 ? 1 : 0;
vector<pair<double,string> > modes =
createWZDecayModes(newTag2,brat,bosons[iother],maxMass);
for(unsigned int ix=0;ix<modes.size();++ix) {
modes[ix].second.replace(modes[ix].second.size() - 1, 1, ";");
createDecayMode(modes[ix].second,modes[ix].first);
brsum += modes[ix].first;
}
}
}
else {
cerr << "BSMModel::readDecay() "
<< "The decay " << tag << " cannot proceed for on-shell\n"
<< "particles, and has too many off-shell gauge bosons,"
<< " skipping it.\n";
}
}
else {
cerr << "BSMModel::readDecay() "
<< "The decay " << tag << " cannot proceed for on-shell\n"
<< "particles, and has too many outgoing gauge bosons skipping it.\n";
}
}
}
}
if( abs(brsum - 1.) > tolerance_ && nmode!=0 ) {
cerr << "Warning: The total branching ratio for " << inpart->PDGName()
<< " from the spectrum file does not sum to 1. The branching fractions"
<< " will be rescaled.\n";
cerr << setprecision(13) << abs(brsum - 1.) << "\n";
}
if(nmode>0) {
inpart->update();
if(inpart->CC()) inpart->CC()->update();
}
}
void BSMModel::createDecayMode(string tag, double brat) const {
tDMPtr dm = generator()->findDecayMode(tag);
if(!dm) {
dm = generator()->preinitCreateDecayMode(tag);
}
generator()->preinitInterface(dm, "OnOff", "set", "On");
generator()->preinitInterface(dm, "Decayer", "set","/Herwig/Decays/Mambo");
ostringstream brf;
brf << setprecision(13)<< brat;
generator()->preinitInterface(dm, "BranchingRatio","set", brf.str());
if(dm->CC()) {
generator()->preinitInterface(dm->CC(), "OnOff", "set", "On");
generator()->preinitInterface(dm->CC(), "Decayer", "set","/Herwig/Decays/Mambo");
generator()->preinitInterface(dm->CC(), "BranchingRatio","set", brf.str());
}
}
vector<pair<double,string> >
BSMModel::createWZDecayModes(string tag, double brat,
tcPDPtr boson, Energy maxMass) const {
vector<pair<double,string> > modes;
double sum(0.);
for(DecaySet::const_iterator dit=boson->decayModes().begin();
dit!=boson->decayModes().end();++dit) {
tcDMPtr mode = *dit;
if(!mode->on()) continue;
string extra;
Energy outMass(ZERO);
for(ParticleMSet::const_iterator pit=mode->products().begin();
pit!=mode->products().end();++pit) {
extra += (**pit).name() + ",";
outMass += (**pit).mass();
}
if(outMass<maxMass) {
sum += mode->brat();
modes.push_back(make_pair(mode->brat(),tag+extra));
}
}
for(unsigned int ix=0;ix<modes.size();++ix)
modes[ix].first *= brat/sum;
return modes;
}
diff --git a/Models/General/BSMModel.h b/Models/General/BSMModel.h
--- a/Models/General/BSMModel.h
+++ b/Models/General/BSMModel.h
@@ -1,162 +1,183 @@
// -*- C++ -*-
#ifndef Herwig_BSMModel_H
#define Herwig_BSMModel_H
//
// This is the declaration of the BSMModel class.
//
#include "Herwig++/Models/StandardModel/StandardModel.h"
#include "ThePEG/Utilities/CFileLineReader.h"
namespace Herwig {
using namespace ThePEG;
/**
* Here is the documentation of the BSMModel class.
*
* @see \ref BSMModelInterfaces "The interfaces"
* defined for BSMModel.
*/
class BSMModel: public Herwig::StandardModel {
public:
/**
* The default constructor.
*/
BSMModel();
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:
/**
* Get name of SLHA decay file
*/
const string & decayFile() const {return decayFile_;}
/**
* Set name of SLHA decay file
*/
void decayFile(string in) {decayFile_ = in;}
/**
* Read the decays
*/
void decayRead();
/**
* Read decaymodes from LHA file
* @param ifs input stream containg data
* @param decay string containing name of parent and value of total width
*/
void readDecay(CFileLineReader & ifs, string decay) const;
/**
* Create a DecayMode object in the repository
* @param tag The tag identifying the decay mode including the prefix
* 'decaymode'
* @param brat Branching ratio of this mode
*/
void createDecayMode(string tag, double brat) const;
/**
* Create a DecayMode object in the repository
* @param tag The tag identifying the decay mode including the prefix
* 'decaymode'
* @param brat Branching ratio of this mode
*/
vector<pair<double,string> > createWZDecayModes(string tag, double brat,
tcPDPtr boson,
Energy maxMass) const;
/**
* read the decays
*/
bool readDecays() const {return readDecays_;}
/**
* set the reading of the decays
*/
void readDecays(bool in) {readDecays_=in;}
+ /**
+ * Map of PDG ids from file to those used internally
+ */
+ map<long,long> & idMap() {return idMap_;}
+
+ /**
+ * Get ParticleData object with Id mapping
+ */
+ PDPtr getBSMParticleData(PID id) const {
+ map<long,long>::const_iterator it = idMap_.find(id);
+ if(it==idMap_.end())
+ return getParticleData(id);
+ else
+ return getParticleData(it->second);
+ }
+
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();
//@}
/**
* Overloaded function from Interfaced
*/
virtual bool preInitialize() const {
return true;
}
private:
/**
* The assignment operator is private and must never be called.
* In fact, it should not even be implemented.
*/
BSMModel & operator=(const BSMModel &);
private:
/**
* Name of the decay file
*/
string decayFile_;
/**
* Read the decays from the file
*/
bool readDecays_;
/**
* Whether or not to replace the top decay modes with those from
* the SLHA files
*/
bool topModesFromFile_;
/**
* Tolerance for branching ratios
*/
double tolerance_;
+ /**
+ * Map of ids from files to those used by Herwig++
+ */
+ map<long,long> idMap_;
+
};
}
#endif /* Herwig_BSMModel_H */
diff --git a/Models/Susy/RPV/RPVFFWVertex.cc b/Models/Susy/RPV/RPVFFWVertex.cc
--- a/Models/Susy/RPV/RPVFFWVertex.cc
+++ b/Models/Susy/RPV/RPVFFWVertex.cc
@@ -1,247 +1,254 @@
// -*- C++ -*-
//
// This is the implementation of the non-inlined, non-templated member
// functions of the RPVFFWVertex class.
//
#include "RPVFFWVertex.h"
#include "ThePEG/Interface/ClassDocumentation.h"
#include "ThePEG/Interface/Switch.h"
#include "ThePEG/EventRecord/Particle.h"
#include "ThePEG/Repository/UseRandom.h"
#include "ThePEG/Repository/EventGenerator.h"
#include "ThePEG/Utilities/DescribeClass.h"
#include "Herwig++/Models/StandardModel/StandardCKM.h"
#include "ThePEG/Persistency/PersistentOStream.h"
#include "ThePEG/Persistency/PersistentIStream.h"
using namespace Herwig;
RPVFFWVertex::RPVFFWVertex() : _ckm(3,vector<Complex>(3,0.0)),
_sw(0.), _couplast(0.), _q2last(ZERO),
_id1last(0), _id2last(0), _leftlast(0.),
_rightlast(0.), _interactions(0) {
orderInGs(0);
orderInGem(1);
}
IBPtr RPVFFWVertex::clone() const {
return new_ptr(*this);
}
IBPtr RPVFFWVertex::fullclone() const {
return new_ptr(*this);
}
void RPVFFWVertex::doinit() {
// SUSY mixing matrices
tSusyBasePtr model = dynamic_ptr_cast<SusyBasePtr>(generator()->standardModel());
if(!model)
throw InitException() << "RPVFFWVertex::doinit() - The model pointer is null!"
<< Exception::abortnow;
_theN = model->neutralinoMix();
_theU = model->charginoUMix();
_theV = model->charginoVMix();
if(!_theN || !_theU || ! _theV)
throw InitException() << "RPVFFWVertex::doinit() - "
<< "A mixing matrix pointer is null."
<< " N: " << _theN << " U: " << _theU << " V: "
<< _theV << Exception::abortnow;
// SM interactions
if(_interactions==0 || _interactions==1) {
// particles for outgoing W-
// quarks
for(int ix=1;ix<6;ix+=2) {
for(int iy=2;iy<7;iy+=2) {
addToList(-ix, iy, -24);
}
}
// leptons
for(int ix=11;ix<17;ix+=2) {
addToList(-ix, ix+1, -24);
}
// particles for outgoing W+
// quarks
for(int ix=2;ix<7;ix+=2) {
for(int iy=1;iy<6;iy+=2) {
addToList(-ix, iy, 24);
}
}
// leptons
for(int ix=11;ix<17;ix+=2) {
addToList(-ix-1, ix, 24);
}
}
// neutralino and chargino
if(_interactions==0 || _interactions==2) {
vector<long> neu(4);
neu[0] = 1000022; neu[1] = 1000023;
neu[2] = 1000025; neu[3] = 1000035;
if(_theN->size().first==7) {
- neu.push_back(12);
- neu.push_back(14);
- neu.push_back(16);
+ if(model->majoranaNeutrinos()) {
+ neu.push_back(17);
+ neu.push_back(18);
+ neu.push_back(19);
+ }
+ else {
+ neu.push_back(12);
+ neu.push_back(14);
+ neu.push_back(16);
+ }
}
vector<long> cha(2);
cha[0] = 1000024; cha[1] = 1000037;
if(_theV->size().first==5) {
cha.push_back(11);
cha.push_back(13);
cha.push_back(15);
}
// sign == -1 outgoing W-, sign == +1 outgoing W+
for(int sign = -1; sign < 2; sign += 2) {
for(unsigned int ine = 0; ine < neu.size(); ++ine) {
for(unsigned int ic = 0; ic < cha.size(); ++ic ) {
if(ic>1&&ine>3&&ic==ine-2) continue;
addToList(-sign*cha[ic], neu[ine], sign*24);
}
}
}
}
Helicity::FFVVertex::doinit();
// CKM matric
Ptr<CKMBase>::transient_pointer CKM = model->CKM();
// cast the CKM object to the HERWIG one
ThePEG::Ptr<Herwig::StandardCKM>::transient_const_pointer
hwCKM = ThePEG::dynamic_ptr_cast< ThePEG::Ptr<Herwig::StandardCKM>::
transient_const_pointer>(CKM);
if(hwCKM) {
vector< vector<Complex > > CKM;
CKM = hwCKM->getUnsquaredMatrix(generator()->standardModel()->families());
for(unsigned int ix=0;ix<3;++ix) {
for(unsigned int iy=0;iy<3;++iy) {
_ckm[ix][iy]=CKM[ix][iy];
}
}
}
else {
throw Exception() << "Must have access to the Herwig::StandardCKM object"
<< "for the CKM matrix in SMFFWVertex::doinit()"
<< Exception::runerror;
}
_sw = sqrt(sin2ThetaW());
}
void RPVFFWVertex::persistentOutput(PersistentOStream & os) const {
os << _sw << _theN << _theU << _theV << _ckm << _interactions;
}
void RPVFFWVertex::persistentInput(PersistentIStream & is, int) {
is >> _sw >> _theN >> _theU >> _theV >> _ckm >> _interactions;
}
// *** Attention *** The following static variable is needed for the type
// description system in ThePEG. Please check that the template arguments
// are correct (the class and its base class), and that the constructor
// arguments are correct (the class name and the name of the dynamically
// loadable library where the class implementation can be found).
DescribeClass<RPVFFWVertex,Helicity::FFVVertex>
describeHerwigRPVFFWVertex("Herwig::RPVFFWVertex", "HwSusy.so HwRPV.so");
void RPVFFWVertex::Init() {
static ClassDocumentation<RPVFFWVertex> documentation
("The couplings of the fermions to the W boson in the RPV model"
" with bilinear R-parity violation");
static Switch<RPVFFWVertex,unsigned int> interfaceInteractions
("Interactions",
"Which interactions to include",
&RPVFFWVertex::_interactions, 0, false, false);
static SwitchOption interfaceInteractionsAll
(interfaceInteractions,
"All",
"Include all the interactions",
0);
static SwitchOption interfaceInteractionsSM
(interfaceInteractions,
"SM",
"Only include the MS terms",
1);
static SwitchOption interfaceInteractionsSUSY
(interfaceInteractions,
"SUSY",
"Include the neutralino/chargino terms",
2);
}
void RPVFFWVertex::setCoupling(Energy2 q2,tcPDPtr part1,
tcPDPtr part2,tcPDPtr part3) {
assert(abs(part3->id()) == ParticleID::Wplus);
// normalization
// first the overall normalisation
if(q2 != _q2last||_couplast==0.) {
_couplast = weakCoupling(q2);
_q2last=q2;
}
norm(_couplast);
// left and right couplings for quarks
if(abs(part1->id()) <= 6) {
int iferm=abs(part1->id());
int ianti=abs(part2->id());
if(iferm%2!=0) swap(iferm,ianti);
iferm = iferm/2;
ianti = (ianti+1)/2;
assert( iferm>=1 && iferm<=3 && ianti>=1 && ianti<=3);
left(-sqrt(0.5)*_ckm[iferm-1][ianti-1]);
right(0.);
}
else {
long neu, cha;
if(part1->charged()) {
cha = part1->id();
neu = part2->id();
}
else {
cha = part2->id();
neu = part1->id();
}
if(_theV->size().first==2&&abs(neu)<=16) {
left(-sqrt(0.5));
right(0.);
}
else {
// assert((abs(cha) == 1000024 || abs(cha) == 1000037) &&
// (neu == 1000022 || neu == 1000023 ||
// neu == 1000025 || neu == 1000035 ||
// neu == 1000045) );
//
if(cha != _id1last || neu != _id2last) {
_id1last = cha;
_id2last = neu;
unsigned int eigc = abs(cha) < 16 ? (abs(cha)-11)/2+2 :
(abs(cha)-1000024)/13;
unsigned int eign = neu<=16 ? 4+(abs(neu)-12)/2 :
(neu<1000025 ? neu - 1000022 : (neu-1000025)/10+2);
_leftlast = (*_theN)(eign, 1)*conj((*_theV)(eigc, 0)) -
( (*_theN)(eign, 3)*conj((*_theV)(eigc, 1))/sqrt(2));
_rightlast = conj((*_theN)(eign, 1))*(*_theU)(eigc, 0) +
( conj((*_theN)(eign, 2))*(*_theU)(eigc, 1)/sqrt(2));
if(_theV->size().first==5) {
for(unsigned int k=0;k<3;++k)
_rightlast += ( conj((*_theN)(eign, 4+k))*(*_theU)(eigc, 2+k)/sqrt(2));
}
}
Complex ltemp = _leftlast;
Complex rtemp = _rightlast;
// conjugate if +ve chargino
if(cha>0) {
ltemp = conj(ltemp);
rtemp = conj(rtemp);
}
if((part1->id()==cha&&cha>0)||(part2->id()==cha&&cha<0)) {
Complex temp = ltemp;
ltemp = -rtemp;
rtemp = -temp;
}
left (ltemp);
right(rtemp);
}
}
}
diff --git a/Models/Susy/RPV/RPVFFZVertex.cc b/Models/Susy/RPV/RPVFFZVertex.cc
--- a/Models/Susy/RPV/RPVFFZVertex.cc
+++ b/Models/Susy/RPV/RPVFFZVertex.cc
@@ -1,277 +1,284 @@
// -*- C++ -*-
//
// This is the implementation of the non-inlined, non-templated member
// functions of the RPVFFZVertex class.
//
#include "RPVFFZVertex.h"
#include "ThePEG/Interface/ClassDocumentation.h"
#include "ThePEG/Interface/Switch.h"
#include "ThePEG/EventRecord/Particle.h"
#include "ThePEG/Repository/UseRandom.h"
#include "ThePEG/Repository/EventGenerator.h"
#include "ThePEG/Utilities/DescribeClass.h"
#include "ThePEG/Persistency/PersistentOStream.h"
#include "ThePEG/Persistency/PersistentIStream.h"
using namespace Herwig;
RPVFFZVertex::RPVFFZVertex() : _sw(0.), _cw(0.), _id1last(0),
_id2last(0), _q2last(), _couplast(0.),
_leftlast(0.), _rightlast(0.),
_gl(17,0.0), _gr(17,0.0), _gblast(0),
_interactions(0) {
orderInGem(1);
orderInGs(0);
}
IBPtr RPVFFZVertex::clone() const {
return new_ptr(*this);
}
IBPtr RPVFFZVertex::fullclone() const {
return new_ptr(*this);
}
void RPVFFZVertex::persistentOutput(PersistentOStream & os) const {
os << _sw << _cw << _theN << _theU << _theV
<< _gl << _gr << _interactions;
}
void RPVFFZVertex::persistentInput(PersistentIStream & is, int) {
is >> _sw >> _cw >> _theN >> _theU >> _theV
>> _gl >> _gr >> _interactions;
}
// *** Attention *** The following static variable is needed for the type
// description system in ThePEG. Please check that the template arguments
// are correct (the class and its base class), and that the constructor
// arguments are correct (the class name and the name of the dynamically
// loadable library where the class implementation can be found).
DescribeClass<RPVFFZVertex,Helicity::FFVVertex>
describeHerwigRPVFFZVertex("Herwig::RPVFFZVertex", "HwSusy.so HwRPV.so");
void RPVFFZVertex::Init() {
static ClassDocumentation<RPVFFZVertex> documentation
("The RPVFFZVertex class implements trhe coupling of the Z to all"
" fermion-antifermion pairs in models with bilinear RPV.");
static Switch<RPVFFZVertex,unsigned int> interfaceInteractions
("Interactions",
"Which interactions to include",
&RPVFFZVertex::_interactions, 0, false, false);
static SwitchOption interfaceInteractionsAll
(interfaceInteractions,
"All",
"Include all the interactions",
0);
static SwitchOption interfaceInteractionsSM
(interfaceInteractions,
"SM",
"Only include what would have been the interactions with the SM"
" fermions in the absence of mixing",
1);
static SwitchOption interfaceInteractionsNeutralino
(interfaceInteractions,
"Neutralino",
"Only include what would have been the interactions with the "
"neutralinos in the absence of mixing",
2);
static SwitchOption interfaceInteractionsChargino
(interfaceInteractions,
"Chargino",
"Only include what would have been the interactions with the "
"charginos in the absence of mixing",
3);
}
void RPVFFZVertex::doinit() {
// extract the mixing matrices
tSusyBasePtr model = dynamic_ptr_cast<SusyBasePtr>(generator()->standardModel());
if(!model) throw InitException() << "RPVFFZVertex::doinit() - "
<< "The model pointer is null."
<< Exception::abortnow;
_theN = model->neutralinoMix();
_theU = model->charginoUMix();
_theV = model->charginoVMix();
if( !_theN || !_theU || !_theV )
throw InitException() << "RPVFFZVertex::doinit - "
<< "A mixing matrix pointer is null. U: "
<< _theU << " V: " << _theV << " N: " << _theN
<< Exception::abortnow;
// Standard Model fermions
if(_interactions==0||_interactions==1) {
// PDG codes for the particles
// the quarks
for(int ix=1;ix<7;++ix) {
addToList(-ix, ix, 23);
}
// the leptons
for(int ix=11;ix<17;ix+=2) {
addToList(-ix, ix, 23);
}
for(int ix=12;ix<17;ix+=2) {
if(_theN->size().first==7)
addToList( ix, ix, 23);
else
addToList(-ix, ix, 23);
}
}
// neutralinos
if(_interactions==0||_interactions==2) {
vector<long> neu(4);
neu[0] = 1000022; neu[1] = 1000023;
neu[2] = 1000025; neu[3] = 1000035;
if(_theN->size().first==7) {
- neu.push_back(12);
- neu.push_back(14);
- neu.push_back(16);
+ if(model->majoranaNeutrinos()) {
+ neu.push_back(17);
+ neu.push_back(18);
+ neu.push_back(19);
+ }
+ else {
+ neu.push_back(12);
+ neu.push_back(14);
+ neu.push_back(16);
+ }
}
for(unsigned int i = 0; i < neu.size(); ++i) {
for(unsigned int j = 0; j < neu.size(); ++j) {
if(i>3&&i==j) addToList(neu[i], neu[j], 23);
}
}
}
// charginos
if(_interactions==0||_interactions==3) {
addToList(-1000024, 1000024, 22);
addToList(-1000037, 1000037, 22);
vector<long> cha(2);
cha[0] = 1000024; cha[1] = 1000037;
if(_theV->size().first==5) {
cha.push_back(11);
cha.push_back(13);
cha.push_back(15);
}
for(unsigned int i = 0; i < cha.size(); ++i) {
for(unsigned int j = 0; j < cha.size(); ++j) {
if(i>1&&i==j) addToList(cha[i], cha[j], 23);
}
}
}
Helicity::FFVVertex::doinit();
// weak mixing
double sw2 = sin2ThetaW();
_cw = sqrt(1. - sw2);
_sw = sqrt( sw2 );
// Standard Model couplings
for(int ix=1;ix<4;++ix) {
_gl[2*ix-1] = -0.25*(model->vd() + model->ad() );
_gl[2*ix ] = -0.25*(model->vu() + model->au() );
_gl[2*ix+9 ] = -0.25*(model->ve() + model->ae() );
_gl[2*ix+10] = -0.25*(model->vnu() + model->anu());
_gr[2*ix-1] = -0.25*(model->vd() - model->ad() );
_gr[2*ix ] = -0.25*(model->vu() - model->au() );
_gr[2*ix+9 ] = -0.25*(model->ve() - model->ae() );
_gr[2*ix+10] = -0.25*(model->vnu() - model->anu());
}
}
void RPVFFZVertex::setCoupling(Energy2 q2,tcPDPtr part1,
tcPDPtr part2,tcPDPtr part3) {
// first the overall normalisation
if(q2!=_q2last||_couplast==0.) {
_couplast = electroMagneticCoupling(q2);
_q2last=q2;
}
long iferm1(part1->id()), iferm2(part2->id()), boson(part3->id());
long iferm = abs(iferm1);
// chargino coupling to the photon
if(part3->id()==ParticleID::gamma) {
assert(iferm == abs(iferm2));
_gblast = boson;
_id1last = iferm1;
_id2last = iferm2;
_leftlast = -1.;
_rightlast = -1.;
}
// coupling to the Z
else {
assert(part3->id()==ParticleID::Z0);
_couplast /= _sw*_cw;
// quarks
if(iferm<=6) {
_leftlast = _gl[iferm];
_rightlast = _gr[iferm];
}
// charged leptons and charginos
else if(part1->iCharge()!=0) {
if(boson != _gblast || iferm1 != _id1last || iferm2 != _id2last) {
_gblast = boson;
_id1last = iferm1;
_id2last = iferm2;
if(_theV->size().first==2&&iferm<=16) {
_leftlast = -_gr[iferm];
_rightlast = -_gl[iferm];
}
else {
unsigned int ic1 = abs(iferm1) < 16 ? (abs(iferm1)-11)/2+2 :
(abs(iferm1)-1000024)/13;
unsigned int ic2 = abs(iferm2) < 16 ? (abs(iferm2)-11)/2+2 :
(abs(iferm2)-1000024)/13;
_leftlast = -(*_theV)(ic1, 0)*conj((*_theV)(ic2, 0)) -
0.5*(*_theV)(ic1, 1)*conj((*_theV)(ic2, 1));
_rightlast = -conj((*_theU)(ic1, 0))*(*_theU)(ic2, 0) -
0.5*conj((*_theU)(ic1, 1))*(*_theU)(ic2, 1);
if(abs(iferm1) == abs(iferm2)) {
_leftlast += sqr(_sw);
_rightlast += sqr(_sw);
}
if(_theV->size().first==5) {
for(unsigned int ix=0;ix<3;++ix) {
_leftlast += -0.5*(*_theV)(ic1, 2+ix)*conj((*_theV)(ic2, 2+ix));
}
}
}
if(iferm1>0) {
Complex temp = _leftlast;
_leftlast = -_rightlast;
_rightlast = -temp;
}
}
}
// neutrinos and neutralinos
else {
// case where only 4x4 matrix and neutrino
if(_theN->size().first==4&&iferm<=16) {
assert(iferm==12||iferm==14||iferm==16);
_leftlast = _gl[iferm];
_rightlast = _gr[iferm];
}
// neutralino
else {
long ic1 = part2->id();
long ic2 = part1->id();
assert(ic1 == ParticleID::SUSY_chi_10 || ic1 == ParticleID::SUSY_chi_20 ||
ic1 == ParticleID::SUSY_chi_30 || ic1 == ParticleID::SUSY_chi_40 ||
abs(ic1) == 12 || abs(ic1) == 14 || abs(ic1) == 16 );
assert(ic2 == ParticleID::SUSY_chi_10 || ic2 == ParticleID::SUSY_chi_20 ||
ic2 == ParticleID::SUSY_chi_30 || ic2 == ParticleID::SUSY_chi_40 ||
abs(ic2) == 12 || abs(ic2) == 14 || abs(ic2) == 16 );
if(ic1 != _id1last || ic2 != _id2last) {
_id1last = ic1;
_id2last = ic2;
unsigned int neu1 = ic1<=16 ? 4+(abs(ic1)-12)/2 :
(ic1<1000025 ? ic1 - 1000022 : (ic1-1000025)/10+2);
unsigned int neu2 = ic2<=16 ? 4+(abs(ic2)-12)/2 :
(ic2<1000025 ? ic2 - 1000022 : (ic2-1000025)/10+2);
_leftlast = 0.5*( (*_theN)(neu1, 3)*conj((*_theN)(neu2, 3)) -
(*_theN)(neu1, 2)*conj((*_theN)(neu2, 2)) );
if(_theN->size().first>4) {
for(unsigned int k=0;k<3;++k)
_leftlast -= 0.5*(*_theN)(neu1, 4+k)*conj((*_theN)(neu2, 4+k));
}
_rightlast = -conj(_leftlast);
}
}
}
}
norm ( _couplast);
left ( _leftlast);
right(_rightlast);
}
diff --git a/Models/Susy/SusyBase.cc b/Models/Susy/SusyBase.cc
--- a/Models/Susy/SusyBase.cc
+++ b/Models/Susy/SusyBase.cc
@@ -1,703 +1,721 @@
// -*- C++ -*-
//
// SusyBase.cc is a part of Herwig++ - A multi-purpose Monte Carlo event generator
// Copyright (C) 2002-2011 The Herwig Collaboration
//
// Herwig++ is licenced under version 2 of the GPL, see COPYING for details.
// Please respect the MCnet academic guidelines, see GUIDELINES for details.
//
//
// This is the implementation of the non-inlined, non-templated member
// functions of the SusyBase class.
//
#include "SusyBase.h"
#include "ThePEG/Utilities/DescribeClass.h"
#include "ThePEG/Interface/ClassDocumentation.h"
#include "ThePEG/Persistency/PersistentOStream.h"
#include "ThePEG/Persistency/PersistentIStream.h"
#include "ThePEG/Interface/Switch.h"
#include "ThePEG/Interface/Reference.h"
#include "ThePEG/Interface/Parameter.h"
#include "ThePEG/Repository/Repository.h"
#include "ThePEG/Utilities/StringUtils.h"
#include "ThePEG/Repository/CurrentGenerator.h"
#include "ThePEG/PDT/EnumParticles.h"
#include "ThePEG/PDT/MassGenerator.h"
#include "ThePEG/PDT/WidthGenerator.h"
#include "ThePEG/PDT/DecayMode.h"
using namespace Herwig;
SusyBase::SusyBase() : readFile_(false), MPlanck_(2.4e18*GeV),
- gravitino_(false),
+ gravitino_(false), majoranaNeutrinos_(false),
tanBeta_(0), mu_(ZERO),
M1_(ZERO), M2_(ZERO), M3_(ZERO),
mH12_(ZERO),mH22_(ZERO),
meL_(ZERO),mmuL_(ZERO),mtauL_(ZERO),
meR_(ZERO),mmuR_(ZERO),mtauR_(ZERO),
mq1L_(ZERO),mq2L_(ZERO),mq3L_(ZERO),
mdR_(ZERO),muR_(ZERO),msR_(ZERO),
mcR_(ZERO),mbR_(ZERO),mtR_(ZERO),
gluinoPhase_(1.)
{}
IBPtr SusyBase::clone() const {
return new_ptr(*this);
}
IBPtr SusyBase::fullclone() const {
return new_ptr(*this);
}
void SusyBase::doinit() {
addVertex(WSFSFVertex_);
addVertex(NFSFVertex_);
addVertex(GFSFVertex_);
addVertex(HSFSFVertex_);
addVertex(CFSFVertex_);
addVertex(GSFSFVertex_);
addVertex(GGSQSQVertex_);
addVertex(GSGSGVertex_);
addVertex(NNZVertex_);
if(NNPVertex_) addVertex(NNPVertex_);
if(GNGVertex_) addVertex(GNGVertex_);
addVertex(CCZVertex_);
addVertex(CNWVertex_);
addVertex(GOGOHVertex_);
addVertex(WHHVertex_);
addVertex(NCTVertex_);
if(gravitino_) {
if(GVNHVertex_) addVertex(GVNHVertex_);
if(GVNVVertex_) addVertex(GVNVVertex_);
if(GVFSVertex_) addVertex(GVFSVertex_);
}
+ if(majoranaNeutrinos()) {
+ idMap().insert(make_pair(12,17));
+ idMap().insert(make_pair(14,18));
+ idMap().insert(make_pair(16,19));
+ }
BSMModel::doinit();
}
void SusyBase::persistentOutput(PersistentOStream & os) const {
os << readFile_ << gravitino_
<< NMix_ << UMix_ << VMix_ << WSFSFVertex_
<< NFSFVertex_ << GFSFVertex_ << HSFSFVertex_ << CFSFVertex_
<< GSFSFVertex_ << GGSQSQVertex_ << GSGSGVertex_
<< NNZVertex_ << NNPVertex_ << CCZVertex_ << CNWVertex_
<< GOGOHVertex_ << WHHVertex_ << GNGVertex_ << NCTVertex_
<< GVNHVertex_ << GVNVVertex_ << GVFSVertex_
<< tanBeta_ << ounit(mu_,GeV)
<< ounit(M1_,GeV) << ounit(M2_,GeV) << ounit(M3_,GeV)
<< ounit(mH12_,GeV2) << ounit(mH22_,GeV2)
<< ounit(meL_,GeV) << ounit(mmuL_,GeV) << ounit(mtauL_,GeV)
<< ounit(meR_,GeV) << ounit(mmuR_,GeV) << ounit(mtauR_,GeV)
<< ounit(mq1L_,GeV) << ounit(mq2L_,GeV) << ounit(mq3L_,GeV)
<< ounit(mdR_,GeV) << ounit(muR_,GeV) << ounit(msR_,GeV)
<< ounit(mcR_,GeV) << ounit(mbR_,GeV) << ounit(mtR_,GeV)
<< gluinoPhase_ << ounit(MPlanck_,GeV);
}
void SusyBase::persistentInput(PersistentIStream & is, int) {
is >> readFile_ >> gravitino_
>> NMix_ >> UMix_ >> VMix_ >> WSFSFVertex_
>> NFSFVertex_ >> GFSFVertex_ >> HSFSFVertex_ >> CFSFVertex_
>> GSFSFVertex_ >> GGSQSQVertex_ >> GSGSGVertex_
>> NNZVertex_ >> NNPVertex_ >> CCZVertex_ >> CNWVertex_
>> GOGOHVertex_ >> WHHVertex_ >> GNGVertex_ >> NCTVertex_
>> GVNHVertex_ >> GVNVVertex_ >> GVFSVertex_
>> tanBeta_ >> iunit(mu_,GeV)
>> iunit(M1_,GeV) >> iunit(M2_,GeV) >> iunit(M3_,GeV)
>> iunit(mH12_,GeV2) >> iunit(mH22_,GeV2)
>> iunit(meL_,GeV) >> iunit(mmuL_,GeV) >> iunit(mtauL_,GeV)
>> iunit(meR_,GeV) >> iunit(mmuR_,GeV) >> iunit(mtauR_,GeV)
>> iunit(mq1L_,GeV) >> iunit(mq2L_,GeV) >> iunit(mq3L_,GeV)
>> iunit(mdR_,GeV) >> iunit(muR_,GeV) >> iunit(msR_,GeV)
>> iunit(mcR_,GeV) >> iunit(mbR_,GeV) >> iunit(mtR_,GeV)
>> gluinoPhase_ >> iunit(MPlanck_,GeV);
}
// Static variable needed for the type description system in ThePEG.
DescribeClass<SusyBase,BSMModel>
describeHerwigSusyBase("Herwig::SusyBase", "HwSusy.so");
void SusyBase::Init() {
static ClassDocumentation<SusyBase> documentation
("This is the base class for any SUSY model.",
"SUSY spectrum files follow the Les Houches accord"
" \\cite{Skands:2003cj,Allanach:2008qq}.",
" %\\cite{Skands:2003cj}\n"
"\\bibitem{Skands:2003cj}\n"
" P.~Skands {\\it et al.},\n"
" ``SUSY Les Houches accord: Interfacing SUSY spectrum calculators, decay\n"
" %packages, and event generators,''\n"
" JHEP {\\bf 0407}, 036 (2004)\n"
" [arXiv:hep-ph/0311123].\n"
" %%CITATION = JHEPA,0407,036;%%\n"
"%\\cite{Allanach:2008qq}\n"
"\\bibitem{Allanach:2008qq}\n"
" B.~Allanach {\\it et al.},\n"
" %``SUSY Les Houches Accord 2,''\n"
" Comput.\\ Phys.\\ Commun.\\ {\\bf 180}, 8 (2009)\n"
" [arXiv:0801.0045 [hep-ph]].\n"
" %%CITATION = CPHCB,180,8;%%\n"
);
-
static Reference<SusyBase,Helicity::AbstractVSSVertex> interfaceVertexWSS
("Vertex/WSFSF",
"Reference to Susy W SF SF vertex",
&SusyBase::WSFSFVertex_, false, false, true, false);
static Reference<SusyBase,Helicity::AbstractFFSVertex> interfaceVertexNFSF
("Vertex/NFSF",
"Reference to the neutralino-fermion-sfermion vertex",
&SusyBase::NFSFVertex_, false, false, true, false);
static Reference<SusyBase,Helicity::AbstractFFSVertex> interfaceVertexGFSF
("Vertex/GFSF",
"Reference to the gluino-fermion-sfermion vertex",
&SusyBase::GFSFVertex_, false, false, true, false);
static Reference<SusyBase,Helicity::AbstractSSSVertex> interfaceVertexHSFSF
("Vertex/HSFSF",
"Reference to the Higgs-fermion-sfermion vertex",
&SusyBase::HSFSFVertex_, false, false, true, false);
static Reference<SusyBase,Helicity::AbstractFFSVertex> interfaceVertexCFSF
("Vertex/CFSF",
"Reference to the chargino-fermion-sfermion vertex",
&SusyBase::CFSFVertex_, false, false, true, false);
static Reference<SusyBase,Helicity::AbstractVSSVertex> interfaceVertexGSFSF
("Vertex/GSFSF",
"Reference to the gluon-sfermion-sfermion vertex",
&SusyBase::GSFSFVertex_, false, false, true, false);
static Reference<SusyBase,Helicity::AbstractVVSSVertex> interfaceVertexGGSS
("Vertex/GGSQSQ",
"Reference to the gluon-gluon-squark-squark vertex",
&SusyBase::GGSQSQVertex_, false, false, true, false);
static Reference<SusyBase,Helicity::AbstractFFVVertex> interfaceVertexGSGSG
("Vertex/GSGSG",
"Reference to the gluon-gluino-gluino vertex",
&SusyBase::GSGSGVertex_, false, false, true, false);
static Reference<SusyBase,Helicity::AbstractFFVVertex> interfaceVertexNNZ
("Vertex/NNZ",
"Reference to Z-~chi_i0-~chi_i0 vertex",
&SusyBase::NNZVertex_, false, false, true, false);
static Reference<SusyBase,Helicity::AbstractFFVVertex> interfaceVertexNNP
("Vertex/NNP",
"Reference to photon-~chi_i0-~chi_i0 vertex",
&SusyBase::NNPVertex_, false, false, true, true, false);
static Reference<SusyBase,Helicity::AbstractFFVVertex> interfaceVertexGNG
("Vertex/GNG",
"Reference to gluon-~chi_i0-gluino vertex",
&SusyBase::GNGVertex_, false, false, true, true, false);
static Reference<SusyBase,Helicity::AbstractFFVVertex> interfaceVertexCCZ
("Vertex/CCZ",
"Reference to ~chi_i+-~chi_i-Z vertex",
&SusyBase::CCZVertex_, false, false, true, false);
static Reference<SusyBase,Helicity::AbstractFFVVertex> interfaceVertexCNW
("Vertex/CNW",
"Reference to ~chi_i+-chi_i0-W vertex",
&SusyBase::CNWVertex_, false, false, true, false);
static Reference<SusyBase,Helicity::AbstractFFSVertex> interfaceVertexGOGOH
("Vertex/GOGOH",
"Reference to the gaugino-gaugino-higgs vertex",
&SusyBase::GOGOHVertex_, false, false, true, false);
static Reference<SusyBase,Helicity::AbstractVSSVertex> interfaceVertexWHH
("Vertex/SSWHH",
"Reference to Susy WHHVertex",
&SusyBase::WHHVertex_, false, false, true, false);
static Reference<SusyBase,AbstractFFSVertex> interfaceVertexNCT
("Vertex/NCT",
"Vertex for the flavour violating coupling of the top squark "
"to the neutralino and charm quark.",
&SusyBase::NCTVertex_, false, false, true, false, false);
static Reference<SusyBase,AbstractRFSVertex> interfaceVertexGVNH
("Vertex/GVNH",
"Vertex for the interfaction of the gravitino-neutralino"
" and Higgs bosons",
&SusyBase::GVNHVertex_, false, false, true, true, false);
static Reference<SusyBase,AbstractRFVVertex> interfaceVertexGVNV
("Vertex/GVNV",
"Vertex for the interfaction of the gravitino-neutralino"
" and vector bosons",
&SusyBase::GVNVVertex_, false, false, true, true, false);
static Reference<SusyBase,AbstractRFSVertex> interfaceVertexGVFS
("Vertex/GVFS",
"Vertex for the interfaction of the gravitino-fermion"
" and sfermion",
&SusyBase::GVFSVertex_, false, false, true, true, false);
static Parameter<SusyBase,Energy> interfaceMPlanck
("MPlanck",
"The Planck mass for GMSB models",
&SusyBase::MPlanck_, GeV, 2.4e18*GeV, 1.e16*GeV, 1.e20*GeV,
false, false, Interface::limited);
+ static Switch<SusyBase,bool> interfaceMajoranaNeutrinos
+ ("MajoranaNeutrinos",
+ "Whether or not the neutrinos should be treated as Majorana particles",
+ &SusyBase::majoranaNeutrinos_, false, false, false);
+ static SwitchOption interfaceMajoranaNeutrinosYes
+ (interfaceMajoranaNeutrinos,
+ "Yes",
+ "Neutrinos are Majorana",
+ true);
+ static SwitchOption interfaceMajoranaNeutrinosNo
+ (interfaceMajoranaNeutrinos,
+ "No",
+ "Neutrinos are Dirac fermions",
+ false);
+
}
void SusyBase::readSetup(istream & is) {
string filename = dynamic_ptr_cast<istringstream*>(&is)->str();
filename = StringUtils::stripws(filename);
if(readFile_)
throw SetupException()
<< "A second SLHA file " << filename << " has been opened."
<< "This is probably unintended and as it can cause crashes"
<< " and other unpredictable behaviour it is not allowed."
<< Exception::runerror;
CFileLineReader cfile;
cfile.open(filename);
if( !cfile ) throw SetupException()
<< "SusyBase::readSetup - An error occurred in opening the "
<< "spectrum file \"" << filename << "\". A SUSY model cannot be "
<< "run without this."
<< Exception::runerror;
useMe();
// read first line and check if this is a Les Houches event file
cfile.readline();
bool lesHouches = cfile.find("<LesHouchesEvents");
bool reading = !lesHouches;
if(lesHouches) cfile.readline();
//function pointer for putting all characters to lower case.
int (*pf)(int) = tolower;
do {
string line = cfile.getline();
// check for start of slha block in SLHA files
if(lesHouches && !reading) {
if(line.find("<slha")==0) reading = true;
if(!cfile.readline()) break;
continue;
}
// ignore comment lines
if(line[0] == '#') {
if(!cfile.readline()) break;
continue;
}
// make everything lower case
transform(line.begin(), line.end(), line.begin(), pf);
// start of a block
if(line.find("block") == 0) {
string name = StringUtils::car(StringUtils::cdr(line), " #");
name = StringUtils::stripws(name);
// mixing matrix
if((name.find("mix") != string::npos &&
name.find("hmix") != 0)) {
unsigned int row(0),col(0);
MixingVector vals = readMatrix(cfile,row,col);
mixings_[name] = make_pair(make_pair(row,col),vals);
}
else if(name.find("au") == 0 || name.find("ad") == 0 ||
name.find("ae") == 0 ) {
string test = StringUtils::car(line, "#");
while (test.find("=")!= string::npos) {
test = StringUtils::cdr(test, "=");
}
istringstream is(test);
double scale;
is >> scale;
unsigned int row(0),col(0);
MixingVector vals = readMatrix(cfile,row,col);
if(scale>1e10) continue;
mixings_[name] = make_pair(make_pair(row,col),vals);
}
else if( name.find("info") == string::npos) {
readBlock(cfile,name,line);
}
else {
if(!cfile.readline()) break;
}
continue;
}
else if( lesHouches && line.find("</slha") == 0 ) {
reading = false;
break;
}
if(!cfile.readline()) break;
}
while(true);
// extract the relevant parameters
extractParameters();
// create the mixing matrices we need
createMixingMatrices();
// set the masses, this has to be done after the
// mixing matrices have been created
resetRepositoryMasses();
// have now read the file
if(decayFile()=="") decayFile(filename);
readFile_=true;
}
void SusyBase::readBlock(CFileLineReader & cfile,string name,string linein) {
if(!cfile)
throw SetupException()
<< "SusyBase::readBlock() - The input stream is in a bad state"
<< Exception::runerror;
// storage or the parameters
string test = StringUtils::car(linein, "#");
ParamMap store;
bool set = true;
// special for the alpha block
if(name.find("alpha") == 0 ) {
double alpha;
cfile.readline();
string line = cfile.getline();
istringstream iss(line);
iss >> alpha;
store.insert(make_pair(1,alpha));
}
else {
// extract the scale from the block if present
if(test.find("=")!= string::npos) {
while(test.find("=")!=string::npos)
test= StringUtils::cdr(test,"=");
istringstream is(test);
double scale;
is >> scale;
// only store the lowest scale block
if(parameters_.find(name)!=parameters_.end()) {
set = scale < parameters_[name][-1];
}
else {
store.insert(make_pair(-1,scale));
}
}
while(cfile.readline()) {
string line = cfile.getline();
// skip comments
if(line[0] == '#') continue;
// reached the end
if( line[0] == 'B' || line[0] == 'b' ||
line[0] == 'D' || line[0] == 'd' ||
line[0] == '<' ) {
cfile.resetline();
break;
}
istringstream is(line);
long index;
double value;
if(name.find("rvlam")!= string::npos||
name=="rvt" || name=="rvtp" || name=="rvtpp") {
int i,j,k;
is >> i >> j >> k >> value;
index = i*100+j*10+k;
}
else {
is >> index >> value;
}
store.insert(make_pair(index, value));
}
}
if(set) parameters_[name]=store;
}
const MixingVector
SusyBase::readMatrix(CFileLineReader & cfile,
unsigned int & row, unsigned int & col) {
if(!cfile)
throw SetupException()
<< "SusyBase::readMatrix() - The input stream is in a bad state."
<< Exception::runerror;
unsigned int rowmax(0), colmax(0);
MixingVector values;
while(cfile.readline()) {
string line = cfile.getline();
// skip comments
if(line[0] == '#') continue;
// reached the end
if( line[0] == 'B' || line[0] == 'b' ||
line[0] == 'D' || line[0] == 'd' ||
line[0] == '<' ) {
cfile.resetline();
break;
}
istringstream is(line);
unsigned int index1, index2;
double real(0.), imag(0.);
is >> index1 >> index2 >> real >> imag;
values.push_back(MixingElement(index1,index2,Complex(real, imag)));
if(index1 > rowmax) rowmax = index1;
if(index2 > colmax) colmax = index2;
}
col=colmax;
row=rowmax;
return values;
}
void SusyBase::createMixingMatrix(MixingMatrixPtr & matrix,
string name, const MixingVector & values,
MatrixSize size) {
matrix = new_ptr(MixingMatrix(size.first,size.second));
for(unsigned int ix=0; ix < values.size(); ++ix)
(*matrix)(values[ix].row-1,values[ix].col-1) = values[ix].value;
vector<long> ids;
if(name == "nmix") {
ids.resize(4);
ids[0] = 1000022; ids[1] = 1000023;
ids[2] = 1000025; ids[3] = 1000035;
}
else if(name == "nmnmix") {
ids.resize(5);
ids[0] = 1000022; ids[1] = 1000023;
ids[2] = 1000025; ids[3] = 1000035;
ids[4] = 1000045;
}
else if(name == "rvnmix") {
ids.resize(7);
ids[0] = 12; ids[1] = 14; ids[2] = 16;
ids[3] = 1000022; ids[4] = 1000023;
- ids[5] = 1000025; ids[6] = 1000035;
+ ids[5] = 1000025; ids[6] = 1000035;
}
else if(name == "rpvmix") {
ids.resize(5);
ids[0] = 1000022; ids[1] = 1000023;
ids[2] = 1000025; ids[3] = 1000035;
ids[4] = 1000045;
}
else if(name == "umix" || name == "vmix") {
ids.resize(2);
- ids[0] = 1000024; ids[1] = 1000037;
+ ids[0] = 1000024; ids[1] = 1000037;
}
else if(name == "rvumix" || name == "rvvmix" ) {
ids.resize(5);
ids[0] = 1000024; ids[1] = 1000037;
ids[2] = -11; ids[3] = -13; ids[4] = -15;
}
else if(name == "stopmix") {
ids.resize(2);
- ids[0] = 1000006; ids[1] = 2000006;
+ ids[0] = 1000006; ids[1] = 2000006;
}
else if(name == "sbotmix") {
ids.resize(2);
- ids[0] = 1000005; ids[1] = 2000005;
+ ids[0] = 1000005; ids[1] = 2000005;
}
else if(name == "staumix") {
ids.resize(2);
- ids[0] = 1000015; ids[1] = 2000015;
+ ids[0] = 1000015; ids[1] = 2000015;
}
else if(name == "nmhmix") {
ids.resize(3);
ids[0] = 25; ids[1] = 35; ids[2] = 45;
}
else if(name == "nmamix") {
ids.resize(2);
ids[0] = 36; ids[1] = 46;
}
else if(name == "rvamix") {
ids.resize(5);
ids[0] = 36; ids[1] = 1000017;
ids[2] = 1000018; ids[3] = 1000019;
ids[4] = 0;
}
else if(name == "rvhmix") {
ids.resize(5);
ids[0] = 25; ids[1] = 35;
- ids[2] = 1000012; ids[3] = 1000014; ids[5] = 1000016;
+ ids[2] = 1000012; ids[3] = 1000014; ids[4] = 1000016;
}
else if(name == "rvlmix") {
ids.resize(8);
ids[0] = 37;
ids[1] = -1000011;
ids[2] = -1000013;
ids[3] = -1000015;
ids[4] = -2000011;
ids[5] = -2000013;
ids[6] = -2000015;
ids[7] = 0;
}
else {
throw SetupException() << "SusyBase::createMixingMatrix() "
<< "cannot find correct title for mixing matrix "
<< name << Exception::runerror;
}
matrix->setIds(ids);
}
void SusyBase::resetRepositoryMasses() {
map<string,ParamMap>::const_iterator fit=parameters_.find("mass");
if(fit==parameters_.end())
throw Exception() << "BLOCK MASS not found in input file"
<< " can't set masses of SUSY particles"
<< Exception::runerror;
ParamMap theMasses = fit->second;
for(ParamMap::iterator it = theMasses.begin(); it != theMasses.end();
++it) {
long id = it->first;
double mass = it->second;
//a negative mass requires an adjustment to the
//associated mixing matrix by a factor of i
if(mass < 0.0) adjustMixingMatrix(id);
PDPtr part = getParticleData(id);
if(!part) throw SetupException()
<< "SusyBase::resetRepositoryMasses() - Particle with PDG code " << id
<< " not found." << Exception::warning;
if(abs(id)<=5||abs(id)==23||abs(id)==24||
(abs(id)>=11&&abs(id)<=16))
cerr << "SusyBase::resetRepositoryMasses() Resetting mass of "
<< part->PDGName() << " using SLHA "
<< "file,\nthis can affect parts of the Standard Model simulation and"
<< " is strongly discouraged.\n";
// reset the masses
resetMass(it->first,it->second*GeV,part);
-
// switch on gravitino interactions?
gravitino_ |= id== ParticleID::SUSY_Gravitino;
}
theMasses.clear();
}
void SusyBase::adjustMixingMatrix(long id) {
//get correct mixing matrix
switch(id) {
case 1000021 :
gluinoPhase_ = Complex(0.,1.);
break;
case 1000022 :
case 1000023 :
case 1000025 :
case 1000035 :
case 1000045 :
case 12 :
case 14 :
case 16 :
if(NMix_)
NMix_->adjustPhase(id);
else
throw SetupException() << "SusyBase::adjustMixingMatrix - "
<< "The neutralino mixing matrix pointer "
<< "is null!" << Exception::runerror;
break;
case 1000024 :
case 1000037 :
if(UMix_)
UMix_->adjustPhase(id);
else
throw SetupException() << "SusyBase::adjustMixingMatrix - "
<< "The U-Type chargino mixing matrix pointer "
<< "is null!" << Exception::runerror;
if(VMix_)
VMix_->adjustPhase(id);
else
throw SetupException() << "SusyBase::adjustMixingMatrix - "
<< "The V-Type chargino mixing matrix pointer "
<< "is null!" << Exception::runerror;
break;
default :
throw SetupException()
<< "SusyBase::adjustMixingMatrix - Trying to adjust mixing matrix "
<< "phase for a particle that does not have a mixing matrix "
<< "associated with it. " << id << " must have a negative mass in "
<< "the spectrum file, this should only occur for particles that mix."
<< Exception::runerror;
}
}
void SusyBase::createMixingMatrices() {
map<string,pair<MatrixSize, MixingVector > >::const_iterator it;
for(it=mixings_.begin();it!=mixings_.end();++it) {
string name=it->first;
// create the gaugino mixing matrices
if(name == "nmix" || name == "nmnmix" || name == "rvnmix" ) {
createMixingMatrix(NMix_,name,it->second.second,it->second.first);
}
else if (name == "umix" || name == "rvumix" ) {
createMixingMatrix(UMix_,name,it->second.second,it->second.first);
}
else if (name == "vmix" || name == "rvvmix" ) {
createMixingMatrix(VMix_,name,it->second.second,it->second.first);
}
}
}
void SusyBase::extractParameters(bool checkmodel) {
map<string,ParamMap>::const_iterator pit;
ParamMap::const_iterator it;
// try and get tan beta from hmin and extpar first
// extract tan beta
tanBeta_ = -1.;
if(tanBeta_<0.) {
pit=parameters_.find("hmix");
if(pit!=parameters_.end()) {
it = pit->second.find(2);
if(it!=pit->second.end()) tanBeta_ = it->second;
}
}
if(tanBeta_<0.) {
pit=parameters_.find("extpar");
if(pit!=parameters_.end()) {
it = pit->second.find(25);
if(it!=pit->second.end()) tanBeta_ = it->second;
}
}
// otherwise from minpar
if(tanBeta_<0.) {
pit=parameters_.find("minpar");
if(pit!=parameters_.end()) {
it = pit->second.find(3);
if(it!=pit->second.end()) tanBeta_ = it->second;
}
}
if(tanBeta_<0.)
throw Exception() << "SusyBase::extractParameters() "
<< "Can't find tan beta in BLOCK MINPAR"
<< " or BLOCK EXTPAR " << Exception::runerror;
if(tanBeta_==0.)
throw Exception() << "Tan(beta) = 0 in SusyBase::extractParameters()"
<< Exception::runerror;
// extract parameters from hmix
pit=parameters_.find("hmix");
if(pit==parameters_.end()) {
if(generator())
generator()->logWarning(Exception("SusyBase::extractParameters() BLOCK HMIX not found setting mu to zero\n",
Exception::warning));
else
cerr << "SusyBase::extractParameters() BLOCK HMIX not found setting mu to zero\n";
mu_=ZERO;
}
else {
mu_=findValue(pit,1,"HMIX","mu")*GeV;
}
pit = parameters_.find("msoft");
if( pit == parameters_.end() )
throw Exception() << "SusyBase::extractParameters() "
<< "BLOCK MSOFT not found in "
<< "SusyBase::extractParameters()"
<< Exception::runerror;
M1_ = findValue(pit,1 ,"MSOFT","M_1" )*GeV;
M2_ = findValue(pit,2 ,"MSOFT","M_2" )*GeV;
M3_ = findValue(pit,3 ,"MSOFT","M_3" )*GeV;
mH12_ = findValue(pit,21,"MSOFT","m_H1^2")*GeV2;
mH22_ = findValue(pit,22,"MSOFT","m_H2^2")*GeV2;
meL_ = findValue(pit,31,"MSOFT","M_eL" )*GeV;
mmuL_ = findValue(pit,32,"MSOFT","M_muL" )*GeV;
mtauL_ = findValue(pit,33,"MSOFT","M_tauL")*GeV;
meR_ = findValue(pit,34,"MSOFT","M_eR" )*GeV;
mmuR_ = findValue(pit,35,"MSOFT","M_muR" )*GeV;
mtauR_ = findValue(pit,36,"MSOFT","M_tauR")*GeV;
mq1L_ = findValue(pit,41,"MSOFT","M_q1L" )*GeV;
mq2L_ = findValue(pit,42,"MSOFT","M_q2L" )*GeV;
mq3L_ = findValue(pit,43,"MSOFT","M_q3L" )*GeV;
muR_ = findValue(pit,44,"MSOFT","M_uR" )*GeV;
mcR_ = findValue(pit,45,"MSOFT","M_cR" )*GeV;
mtR_ = findValue(pit,46,"MSOFT","M_tR" )*GeV;
mdR_ = findValue(pit,47,"MSOFT","M_dR" )*GeV;
msR_ = findValue(pit,48,"MSOFT","M_sR" )*GeV;
mbR_ = findValue(pit,49,"MSOFT","M_bR" )*GeV;
if(checkmodel) {
throw Exception() << "The SusyBase class should not be used as a "
<< "Model class, use one of the models which inherit"
<< " from it" << Exception::runerror;
}
}
diff --git a/Models/Susy/SusyBase.h b/Models/Susy/SusyBase.h
--- a/Models/Susy/SusyBase.h
+++ b/Models/Susy/SusyBase.h
@@ -1,669 +1,679 @@
// -*- C++ -*-
//
// SusyBase.h is a part of Herwig++ - A multi-purpose Monte Carlo event generator
// Copyright (C) 2002-2011 The Herwig Collaboration
//
// Herwig++ is licenced under version 2 of the GPL, see COPYING for details.
// Please respect the MCnet academic guidelines, see GUIDELINES for details.
//
#ifndef HERWIG_SusyBase_H
#define HERWIG_SusyBase_H
//
// This is the declaration of the SusyBase class.
//
#include "Herwig++/Models/General/BSMModel.h"
#include "MixingMatrix.h"
#include "ThePEG/Utilities/CFileLineReader.h"
#include "ThePEG/Helicity/Vertex/AbstractVSSVertex.h"
#include "ThePEG/Helicity/Vertex/AbstractSSSVertex.h"
#include "ThePEG/Helicity/Vertex/AbstractVVSSVertex.h"
#include "ThePEG/Helicity/Vertex/AbstractRFSVertex.h"
#include "ThePEG/Helicity/Vertex/AbstractRFVVertex.h"
#include "SusyBase.fh"
namespace Herwig {
using namespace ThePEG;
/*@name Some convenient typedefs. */
//@{
/**
* Map to hold key, parameter pairs.
*/
typedef map<long, double> ParamMap;
//@}
/** \ingroup Models
* This class is designed to be a base class for SUSY models. There is
* an interface to set the name of the spectrum file to read in
* the necessary parameters for a SUSY model.
*
* @see \ref SusyBaseInterfaces "The interfaces"
* defined for SusyBase.
* @see StandardModel
*/
class SusyBase: public BSMModel {
public:
/**
* The default constructor.
*/
SusyBase();
public:
/** @name Access functions. */
//@{
/**
* Value of \f$\tan\beta\f$.
*/
double tanBeta() const { return tanBeta_; }
/**
* Value of \f$\mu\f$ parameter.
*/
Energy muParameter() const { return mu_; }
/**
* The neutralino mixing matrix
*/
const MixingMatrixPtr & neutralinoMix() const {
return NMix_;
}
/**
* The U-type chargino mixing matrix
*/
const MixingMatrixPtr & charginoUMix() const {
return UMix_;
}
/**
* The V-type chargino mixing matrix
*/
const MixingMatrixPtr & charginoVMix() const {
return VMix_;
}
/**
* The phase for gluino vertices
*/
const Complex & gluinoPhase() const {return gluinoPhase_;}
//@}
+ /**
+ * Treatment of neutrinos
+ */
+ bool majoranaNeutrinos() const {return majoranaNeutrinos_;}
+
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();
public:
/**
* Soft breaking parameters
*/
//@{
/**
* The bilinear breaking mass term for the bino
*/
const Energy & M1() const {return M1_;}
/**
* The bilinear breaking mass term for the wino
*/
const Energy & M2() const {return M2_;}
/**
* The bilinear breaking mass term for the gluinos
*/
const Energy & M3() const {return M3_;}
/**
* The soft breaking mass squared for \f$H_1\f$
*/
const Energy2 & Mh12() const {return mH12_;}
/**
* The soft breaking mass squared for \f$H_2\f$
*/
const Energy2 & Mh22() const {return mH22_;}
/**
* Soft breaking mass for the first generation lepton doublet
*/
const Energy & MeL() const {return meL_;}
/**
* Soft breaking mass for the second generation lepton doublet
*/
const Energy & MmuL() const {return mmuL_;}
/**
* Soft breaking mass for the third generation lepton doublet
*/
const Energy & MtauL() const {return mtauL_;}
/**
* Soft breaking mass for the first generation lepton singlet
*/
const Energy & MeR() const {return meR_;}
/**
* Soft breaking mass for the second generation lepton singlet
*/
const Energy & MmuR() const {return mmuR_;}
/**
* Soft breaking mass for the third generation lepton singlet
*/
const Energy & MtauR() const {return mtauR_;}
/**
* Soft breaking mass for the first generation quark doublet
*/
const Energy & Mq1L() const {return mq1L_;}
/**
* Soft breaking mass for the second generation quark doublet
*/
const Energy & Mq2L() const {return mq2L_;}
/**
* Soft breaking mass for the third generation quark doublet
*/
const Energy & Mq3L() const {return mq3L_;}
/**
* Soft breaking mass for the down singlet
*/
const Energy & MdR() const {return mdR_;}
/**
* Soft breaking mass for the up singlet
*/
const Energy & MuR() const {return muR_;}
/**
* Soft breaking mass for the strange singlet
*/
const Energy & MsR() const {return msR_;}
/**
* Soft breaking mass for the charm singlet
*/
const Energy & McR() const {return mcR_;}
/**
* Soft breaking mass for the bottom singlet
*/
const Energy & MbR() const {return mbR_;}
/**
* Soft breaking mass for the top singlet
*/
const Energy & MtR() const {return mtR_;}
//@}
/**
* Planck mass
*/
const Energy & MPlanck() const {return MPlanck_;}
protected:
/**
* Function to read information from a setup file.
* @param is istream object to read file.
*/
virtual void readSetup(istream & is);
private:
/**@name Functions to help file read-in. */
//@{
/**
* Read block from LHA file
* @param ifs input stream containg data
* @param name The name of the block
* @param line The line defining the block
*/
void readBlock(CFileLineReader & ifs,string name,string line);
/**
* Function to read mixing matrix from LHA file
* @param ifs input stream containg data
* @param row Number of rows
* @param col Number of columns
*/
const MixingVector readMatrix(CFileLineReader & ifs, unsigned int & row,
unsigned int & col);
protected:
/**
* Create the mixing matrices for the model
*/
virtual void createMixingMatrices();
/**
* Extract the parameters from the input blocks
*/
virtual void extractParameters(bool checkmodel=true);
/**
* Create a object MixingMatrix in the repository
* @param matrix Pointer to the mixing matrix
* @param name Name of the mixing matrix, i.e. nmix, umix...
* @param values Value of each entry in the matrix
* @param size The size of the matrix
*/
void createMixingMatrix(MixingMatrixPtr & matrix, string name,
const MixingVector & values,
MatrixSize size);
/**
* Reset masses in the repository to values read from LHA file.
*/
void resetRepositoryMasses();
/**
* Adjust row of Mixing Matrix if a negative mass occurs in LHA file
* @param id The PDG code of the particle with a negative mass
*/
virtual void adjustMixingMatrix(long id);
//@}
/**
* Access to the mixings and parameters for the inheriting classes
*/
//@{
/**
* Parameter blocks
*/
const map<string,ParamMap> & parameters() const {
return parameters_;
}
/**
* Mixing blocks
*/
const map<string,pair<MatrixSize,MixingVector> > & mixings() const {
return mixings_;
}
//@}
/**
* Reset neutralino mixing matrix
*/
void neutralinoMix(MixingMatrixPtr nm) { NMix_ = nm; }
/**
* Reset the U-type chargino mixing matrix
*/
void charginoUMix(MixingMatrixPtr um) { UMix_ = um; }
/**
* Reset the V-type chargino mixing matrix
*/
void charginoVMix(MixingMatrixPtr vm) { VMix_ = vm; }
/**
* Read a parameter from a block, checking that the
* entry exists
*/
double findValue(const map<string,ParamMap>::const_iterator pit,
int iloc, const string & block,
const string & name) {
ParamMap::const_iterator it = pit->second.find(iloc);
if(it!=pit->second.end()) {
return it->second;
}
else {
ostringstream message;
message << "SusyBase::findValue() Parameter " << name << " = " << iloc
<< " not found in BLOCK " << block << "\n";
if(generator())
generator()->logWarning( Exception(message.str(), Exception::warning) );
else
cerr << message.str();
return 0.;
}
}
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:
/**
* The assignment operator is private and must never be called.
* In fact, it should not even be implemented.
*/
SusyBase & operator=(const SusyBase &);
private:
/**
- * Whether or not the SLHA fiel has been read
+ * Whether or not the SLHA file has been read
*/
bool readFile_;
/**
* Planck mass needed in GMSB models
*/
Energy MPlanck_;
/**
* Whether or not to include gravitino interactions
*/
bool gravitino_;
+ /**
+ * Treatment of the neutrinos
+ */
+ bool majoranaNeutrinos_;
+
/*
* Storage of the parameters.
*/
//@{
/**
* Parameter blocks
*/
map<string,ParamMap> parameters_;
/**
* Mixing blocks
*/
map<string,pair<MatrixSize, MixingVector> > mixings_;
/**
* \f$\tan\beta\f$
*/
double tanBeta_;
/**
* \f$\mu\f$
*/
Energy mu_;
//@}
/**
* Soft breaking parameters
*/
//@{
/**
* The bilinear breaking mass term for the bino
*/
Energy M1_;
/**
* The bilinear breaking mass term for the wino
*/
Energy M2_;
/**
* The bilinear breaking mass term for the gluinos
*/
Energy M3_;
/**
* The soft breaking mass squared for \f$H_1\f$
*/
Energy2 mH12_;
/**
* The soft breaking mass squared for \f$H_2\f$
*/
Energy2 mH22_;
/**
* Soft breaking mass for the first generation lepton doublet
*/
Energy meL_;
/**
* Soft breaking mass for the second generation lepton doublet
*/
Energy mmuL_;
/**
* Soft breaking mass for the third generation lepton doublet
*/
Energy mtauL_;
/**
* Soft breaking mass for the first generation lepton singlet
*/
Energy meR_;
/**
* Soft breaking mass for the second generation lepton singlet
*/
Energy mmuR_;
/**
* Soft breaking mass for the third generation lepton singlet
*/
Energy mtauR_;
/**
* Soft breaking mass for the first generation quark doublet
*/
Energy mq1L_;
/**
* Soft breaking mass for the second generation quark doublet
*/
Energy mq2L_;
/**
* Soft breaking mass for the third generation quark doublet
*/
Energy mq3L_;
/**
* Soft breaking mass for the down singlet
*/
Energy mdR_;
/**
* Soft breaking mass for the up singlet
*/
Energy muR_;
/**
* Soft breaking mass for the strange singlet
*/
Energy msR_;
/**
* Soft breaking mass for the charm singlet
*/
Energy mcR_;
/**
* Soft breaking mass for the bottom singlet
*/
Energy mbR_;
/**
* Soft breaking mass for the top singlet
*/
Energy mtR_;
//@}
/**
* Phase for the gluino
*/
Complex gluinoPhase_;
/**
* Neutralino and Chargino mixing matrices
*/
//@{
/**
* Store pointers to the gaugino mixing matrices
*/
//@{
/**
* The neutralino mixing matrix
*/
MixingMatrixPtr NMix_;
/**
* The \f$U\f$ mixing matrix for the charginos
*/
MixingMatrixPtr UMix_;
/**
* The \f$V\f$ mixing matrix for the charginos
*/
MixingMatrixPtr VMix_;
//@}
/**@name Vertex pointers. */
//@{
/**
* Pointer to the gauge boson sfermion-sfermion vertex
*/
AbstractVSSVertexPtr WSFSFVertex_;
/**
* Pointer to the neutralino-fermion-sfermion vertex
*/
AbstractFFSVertexPtr NFSFVertex_;
/**
* Pointer to the gluino-fermion-sfermion coupling
*/
AbstractFFSVertexPtr GFSFVertex_;
/**
* Pointer to the Higgs-sfermion-sfermion vertex
*/
AbstractSSSVertexPtr HSFSFVertex_;
/**
* Pointer to the \f$\tilde{\chi}^+\f$-fermion-sfermion vertex
*/
AbstractFFSVertexPtr CFSFVertex_;
/**
* Pointer to the gluon-sfermion-sfermion vertex
*/
AbstractVSSVertexPtr GSFSFVertex_;
/**
* Pointer to the gluon-gluon-squark-squark vertex;
*/
AbstractVVSSVertexPtr GGSQSQVertex_;
/**
* Pointer to the gluon-gluino-gluino vertex
*/
AbstractFFVVertexPtr GSGSGVertex_;
/**
* Pointer to the gluino-neutralino-gluon vertex
*/
AbstractFFVVertexPtr GNGVertex_;
/**
* Pointer to the neutralino-neutralino-Z vertex
*/
AbstractFFVVertexPtr NNZVertex_;
/**
* Pointer to the neutralino-neutralino-photon vertex
*/
AbstractFFVVertexPtr NNPVertex_;
/**
* Pointer to the vertex chargino-chargino-Z vertex
*/
AbstractFFVVertexPtr CCZVertex_;
/**
* Pointer to the vertex chargino-neutralino-Z vertex
*/
AbstractFFVVertexPtr CNWVertex_;
/**
* Pointer to the vertex gaugino-gaugino-higgs vertex
*/
AbstractFFSVertexPtr GOGOHVertex_;
/**
* Pointer to the vertex for a gauge boson and higgs
*/
AbstractVSSVertexPtr WHHVertex_;
/**
* Pointer to the vertex for flavour changing stop decay
*/
AbstractFFSVertexPtr NCTVertex_;
/**
* Pointer to the vertex for gravitino-neutralino vector boson
*/
AbstractRFVVertexPtr GVNVVertex_;
/**
* Pointer to the vertex for gravitino-neutralino Higgs boson
*/
AbstractRFSVertexPtr GVNHVertex_;
/**
* Pointer to the vertex for gravitino-fermion sfermion
*/
AbstractRFSVertexPtr GVFSVertex_;
//@}
};
}
#endif /* HERWIG_SusyBase_H */
diff --git a/src/RPV.model b/src/RPV.model
--- a/src/RPV.model
+++ b/src/RPV.model
@@ -1,92 +1,101 @@
##################################################
# Common setup for the RPV MSSM
#
###################################################
#
# Create particle content
#
###################################################
# load MSSM first
read MSSM.model
+# Create Majorana neutrinos
+cd /Herwig/Particles
+create ThePEG::ParticleData nu_1
+setup nu_1 17 nu_1 0 0. 0 0 0 0 2 1
+create ThePEG::ParticleData nu_2
+setup nu_2 18 nu_2 0 0. 0 0 0 0 2 1
+create ThePEG::ParticleData nu_3
+setup nu_3 19 nu_3 0 0. 0 0 0 0 2 1
+
###################################################
#
# Main directory and model object
#
###################################################
mkdir /Herwig/NewPhysics/RPV
cd /Herwig/NewPhysics/RPV
create Herwig::RPV Model HwRPV.so
# SM couplings
set Model:QCD/RunningAlphaS /Herwig/AlphaS
set Model:EW/RunningAlphaEM /Herwig/AlphaEM
set Model:EW/CKM /Herwig/CKM
set Model:RunningMass /Herwig/RunningMass
###################################################
#
# Vertices
#
###################################################
# create RPV MSSM vertices
mkdir /Herwig/Vertices/RPV
cd /Herwig/Vertices/RPV
create Herwig::RPVLLEVertex RPV_LLE
create Herwig::RPVLQDVertex RPV_LQD
create Herwig::RPVUDDVertex RPV_UDD
# set the vertices
cd /Herwig/NewPhysics/RPV
# SM vertices
set Model:Vertex/FFZ /Herwig/Vertices/FFZVertex
set Model:Vertex/FFW /Herwig/Vertices/FFWVertex
set Model:Vertex/FFG /Herwig/Vertices/FFGVertex
set Model:Vertex/FFP /Herwig/Vertices/FFPVertex
set Model:Vertex/GGG /Herwig/Vertices/GGGVertex
set Model:Vertex/GGGG /Herwig/Vertices/GGGGVertex
set Model:Vertex/WWW /Herwig/Vertices/WWWVertex
set Model:Vertex/WWWW /Herwig/Vertices/WWWWVertex
# MSSM feynman rules
set Model:Vertex/WSFSF /Herwig/Vertices/MSSM/MSSM_WSS
set Model:Vertex/NFSF /Herwig/Vertices/MSSM/MSSM_NFS
set Model:Vertex/GFSF /Herwig/Vertices/MSSM/MSSM_GFS
set Model:Vertex/HSFSF /Herwig/Vertices/MSSM/MSSM_HSS
set Model:Vertex/CFSF /Herwig/Vertices/MSSM/MSSM_CFS
set Model:Vertex/GSFSF /Herwig/Vertices/MSSM/MSSM_GSS
set Model:Vertex/GGSQSQ /Herwig/Vertices/MSSM/MSSM_GGSS
set Model:Vertex/GSGSG /Herwig/Vertices/MSSM/MSSM_GGOGO
set Model:Vertex/GNG /Herwig/Vertices/MSSM/MSSM_GNG
set Model:Vertex/NNZ /Herwig/Vertices/MSSM/MSSM_NNZ
set Model:Vertex/NNP /Herwig/Vertices/MSSM/MSSM_NNP
set Model:Vertex/CCZ /Herwig/Vertices/MSSM/MSSM_CCZ
set Model:Vertex/CNW /Herwig/Vertices/MSSM/MSSM_CNW
set Model:Vertex/FFH /Herwig/Vertices/MSSM/MSSM_FFH
set Model:Vertex/GOGOH /Herwig/Vertices/MSSM/MSSM_GOGOH
set Model:Vertex/WWH /Herwig/Vertices/MSSM/MSSM_WWH
set Model:Vertex/SSWHH /Herwig/Vertices/MSSM/MSSM_WHH
set Model:Vertex/HHH /Herwig/Vertices/MSSM/MSSM_HHH
set Model:Vertex/HGG /Herwig/Vertices/MSSM/MSSM_HGG
set Model:Vertex/HPP /Herwig/Vertices/MSSM/MSSM_HPP
set Model:Vertex/WWHH /Herwig/Vertices/MSSM/MSSM_WWHH
set Model:Vertex/NCT /Herwig/Vertices/MSSM/MSSM_NCT
set Model:Vertex/GVNH /Herwig/Vertices/MSSM/MSSM_GVNH
set Model:Vertex/GVFS /Herwig/Vertices/MSSM/MSSM_GVFS
set Model:Vertex/GVNV /Herwig/Vertices/MSSM/MSSM_GVNV
# RPV Feynman rules
set Model:Vertex/LLE /Herwig/Vertices/RPV/RPV_LLE
set Model:Vertex/LQD /Herwig/Vertices/RPV/RPV_LQD
set Model:Vertex/UDD /Herwig/Vertices/RPV/RPV_UDD
###################################################
# Set up the bsm framework
###################################################
cd /Herwig/NewPhysics
set RPV/Model:ModelGenerator NewModel
set TwoBodyDC:CreateDecayModes Yes
set ThreeBodyDC:CreateDecayModes Yes
###################################################
#
# Choose RPV MSSM over SM
#
###################################################
cd /Herwig/Generators
set LEPGenerator:StandardModelParameters /Herwig/NewPhysics/RPV/Model
set LHCGenerator:StandardModelParameters /Herwig/NewPhysics/RPV/Model

File Metadata

Mime Type
text/x-diff
Expires
Mon, Jan 20, 10:17 PM (1 d, 10 h)
Storage Engine
blob
Storage Format
Raw Data
Storage Handle
4215746
Default Alt Text
(77 KB)

Event Timeline