Page MenuHomeHEPForge
Files F7877257

No OneTemporary
Actions

View Options

diff --git a/mcfm/mcfm_interface.C b/mcfm/mcfm_interface.C
index 60e478b..e2bcc5d 100644
--- a/mcfm/mcfm_interface.C
+++ b/mcfm/mcfm_interface.C
@@ -1,193 +1,191 @@
#include "mcfm_interface.h"
#include "settings.h"
#include "coupling.h"
#include "resconst.h"
#include "phasespace.h"
#include <cstring>
#include <iostream>
#include <math.h>
map <int,string> plabel;
void mcfm::init()
{
energy_.sroot_ = opts.sroot;
density_.ih1_ = opts.ih1;
density_.ih2_ = opts.ih2;
nproc_.nproc_ = opts.nproc;
zerowidth_.zerowidth_ = opts.zerowidth;
dynamicscale_.dynamicscale_ = opts.dynamicscale; // Dynamic scale (if true mufac=muren=q)
// CKM matrix entries
cabib_.Vud_ = opts.Vud;
cabib_.Vus_ = opts.Vus;
cabib_.Vub_ = opts.Vub;
cabib_.Vcd_ = opts.Vcd;
cabib_.Vcs_ = opts.Vcs;
cabib_.Vcb_ = opts.Vcb;
dymasses_.wwidth_ = opts.wwidth;
dymasses_.zwidth_ = opts.zwidth;
//initialise MCFM settings
- flag_.flag_ = false;
-
noglue_.noglue_=false;
noglue_.ggonly_=false;
noglue_.gqonly_=false;
string part = "tota";
strncpy(part_.part_ , part.c_str(), part.size());
colc_.colourchoice_=0;
double rtsmin = 40;
cutoff_.cutoff_=0.001; //add to settings (this is a cut off in mcfm on invariant mass pairs between emitted and radiator)
flags_.qflag_=true;
flags_.gflag_=true;
//Catani-Seymour subtraction cut-offs for initial-initial, initial-final, final-initial, and final-final dipoles
//(add to settings)
alfacut_.aii_=1.;//0.01;
alfacut_.aif_=1.;//0.01;
alfacut_.afi_=1.;//0.01;
alfacut_.aff_=1.;//0.01;
//Limits on invariant mass of vector boson decay products
//(irrelevant if zerowidth=true)
limits_.wsqmin_=pow(phasespace::mmin,2);
limits_.wsqmax_=pow(phasespace::mmax,2);
//Check if the limits are compatible with sroot
if (limits_.wsqmax_> pow(energy_.sroot_,2))
limits_.wsqmax_=pow(energy_.sroot_,2);
plabel[1]="pp";
plabel[2]="pp";
//dycoupling_();
coupling::init();
dymasses_.mb_=0; //probably irrelevant to set the bottom mass to 0, since it is not used anywhere
if (!((density_.ih1_==1) && (density_.ih2_ == -1))
&& !((density_.ih1_==1) && (density_.ih2_ == 1)))
{
cout << "The selected initial state is not valid," << endl;
cout << "please select proton-proton (1 1) or proton-antiproton (1 -1)" << endl;
exit(0);
}
//the default behaviour is to remove no branching ratio
bool removebr=false;
brnrat_.brnrat_=1.;
//branching ratios
double brwen,brzee,brtau,brtop;
branch_(brwen,brzee,brtau,brtop);
//W+->l+nubar
if (nproc_.nproc_ == 1)
{
plabel[3]="nl";
plabel[4]="ea";
plabel[5]="pp";
plabel[6]="pp";
nwz_.nwz_=1;
breit_.mass3_= dymasses_.wmass_;
breit_.width3_= dymasses_.wwidth_;
opts.rmass = dymasses_.wmass_;
opts.rwidth = dymasses_.wwidth_;
if (removebr)
brnrat_.brnrat_=brwen;
}
//W-=>l-nu
else if(nproc_.nproc_ == 2)
{
plabel[3]="el";
plabel[4]="na";
plabel[5]="pp";
plabel[6]="pp";
nwz_.nwz_=-1;
breit_.mass3_=dymasses_.wmass_;
breit_.width3_=dymasses_.wwidth_;
opts.rmass = dymasses_.wmass_;
opts.rwidth = dymasses_.wwidth_;
if (removebr)
brnrat_.brnrat_=brwen;
}
else if(nproc_.nproc_==3)
{
plabel[3]="el";
plabel[4]="ea";
plabel[5]="pp";
plabel[6]="pp";
nwz_.nwz_=0;
breit_.mass3_=dymasses_.zmass_;
breit_.width3_=dymasses_.zwidth_;
opts.rmass = dymasses_.zmass_;
opts.rwidth = dymasses_.zwidth_;
zcouple_.l1_=zcouple_.le_;
zcouple_.r1_=zcouple_.re_;
//with q1=0 the photon contribution is switched off
//with q1=-1 the photon contribution is switched on
dyphoton_.phot_=zcouple_.q1_;
if (removebr)
brnrat_.brnrat_=brzee;
}
else
{
cout << "Wrong process" << endl;
exit(0);
}
//Breit Wigner unweighting settings
//the jets, particles 5 and 6, are not resonant
breit_.n2_=0;
breit_.mass2_=0;
breit_.width2_=0;
//the leptons, particles 3 and 4, are resonant
breit_.n3_=1;
//Decide if worth using Breit-Wigner unweighting or not for the dilepon system
//n3=0
//vmass=wmass
//if(nproc.eq.3) vmass=zmass
//if(mwmin.lt.vmass.and.mwmax.gt.vmass) n3=1
ckmfill_(nwz_.nwz_);
// Set-up PS integration cut-offs
rtsmin = min (rtsmin, sqrt(limits_.wsqmin_ + cutoff_.cutoff_));
if (zerowidth_.zerowidth_)
{
rtsmin = dymasses_.wmass_;
if (nproc_.nproc_ == 3) rtsmin = dymasses_.zmass_;
}
taumin_.taumin_=pow((rtsmin/energy_.sroot_),2);
taumin_.logtaumin_=log(taumin_.taumin_);
xmin_.xmin_=taumin_.taumin_;
//Set-up incoming beams (used only in realint.f)
pext_.p1ext_[3]=-0.5*energy_.sroot_;
pext_.p1ext_[0]=0.;
pext_.p1ext_[1]=0.;
pext_.p1ext_[2]=-0.5*energy_.sroot_;
pext_.p2ext_[3]=-0.5*energy_.sroot_;
pext_.p2ext_[0]=0.;
pext_.p2ext_[1]=0.;
pext_.p2ext_[2]=+0.5*energy_.sroot_;
//set NF to 5 (it is used in H2calc)
nf_.nf_ = resconst::NF;
}
diff --git a/src/init.C b/src/init.C
index e21b3ee..94a8f96 100644
--- a/src/init.C
+++ b/src/init.C
@@ -1,84 +1,85 @@
#include "banner.h"
#include "interface.h"
#include "mcfm_interface.h"
#include "settings.h"
#include "init.h"
#include "pdf.h"
#include "pdfevol.h"
#include "pegasus.h"
#include "coupling.h"
#include "gaussrules.h"
#include "mellinint.h"
#include "rapint.h"
#include "mesq.h"
#include "resconst.h"
#include "anomalous.h"
#include "resint.h"
#include "vjint.h"
#include "vjloint.h"
#include "switch.h"
#include "plotter.h"
#include "printsettings.h"
#include "cubacall.h"
#include <cuba.h>
#include <math.h>
#include <iostream>
#include <cstring>
//rewritten initialisation functions
void dyturboinit(int argc, char * argv[])
{
banner();
gaussinit_(); //initialisation of fortran gaussian quadrature nodes and weights
coupling::SMparameters(); //initialisation of unused MCFM parameters
opts.parse_options(argc,argv);
dofill_.doFill_ = 0;
//Initialise some DYRES settings
g_param_.g_param_ = opts.g_param;
nnlo_.order_ = opts.order; //order (0=LO, 1=NLO, 2=NNLO)
opts_.fixedorder_ = opts.fixedorder; //fixed order/resummation switch
qtsub_.xqtcut_= opts.xqtcut; //Cut on qt/Q
qtsub_.qtcut_= opts.qtcut; //Cut on qt
//move here the flaq.eq.0 initialisation part of resumm() in main2 instead of using this initialisation flag
+ flag_.flag_ = false;
mcfm::init();
iniflavreduce_(); //need to call this after nproc_.nproc_ is set
coupling::initscales();
//C++ resum
//initialise all the C modules
gr::init(); //nodes and weights of gaussian quadrature rules
mellinint::initgauss(); //gaussian quadrature for mellin inversion
mesq::init(); //EW couplings for born amplitudes
rapint::init(); //allocate memory for the rapidity quadrature
resconst::init(); //calculate beta, A and B coefficients
anomalous::init(); //calculate anomalous dimensions, C1, C2 and gamma coefficients
pdfevol::init(); //transform the PDF from x- to N-space at the factorisation scale
pegasus::init(); //initialise Pegasus QCD and transform the PDF from x- to N-space at the starting scale
resint::init(); //initialise dequad integration for the bessel integral
//end C++ resum
//V+j fixed order initialisation
vjint::init();
vjloint::init();
switching::init(); //switching function initialisation
rescinit_();
cubacores(opts.cubacores,1000000); //< set number of cores (move this to cubainit)
cubaexit((void (*)()) exitfun,NULL); //< merge at the end of the run
// histogram output
hists.Init();
/***********************************/
//print out EW and QCD parameters and other settings
if (opts.verbose)
opts.dumpAll();
printsettings();
/***********************************/
}

File Metadata

Mime Type
text/x-diff
Expires
Tue, Nov 19, 2:58 PM (1 d, 16 h)
Storage Engine
blob
Storage Format
Raw Data
Storage Handle
3804886
Default Alt Text
(8 KB)

Event Timeline