diff --git a/GeneratorInterface/Pythia8Interface/plugins/PowhegHooksBB4L.h b/GeneratorInterface/Pythia8Interface/plugins/PowhegHooksBB4L.h index 411273ed0d612..781d36c87209a 100644 --- a/GeneratorInterface/Pythia8Interface/plugins/PowhegHooksBB4L.h +++ b/GeneratorInterface/Pythia8Interface/plugins/PowhegHooksBB4L.h @@ -1,148 +1,195 @@ // PowhegHooksBB4L.h -// Copyright (C) 2017 Silvia Ferrario Ravasio, Tomas Jezo, Paolo Nason, Markus Seidel -// inspired by PowhegHooks.h by Richard Corke -// adjusted to work with EmissionVetoHook1 in CMSSW by Alexander Grohsjean +// Rewritten by T. Jezo in 2021. With various contributions from S. Ferrario +// Ravasio, B. Nachman, P. Nason and M. Seidel. Inspired by +// ttb_NLO_dec/main-PYTHIA8.f by P. Nason and E. Re and by PowhegHooks.h by R. +// Corke. +// +// Adapted for CMSSW by Laurids Jeppe. + +// # Introduction +// +// This hook is intended for use together with POWHEG-BOX-RES/b_bbar_4l NLO LHE +// events. This also includes events in which one of the W bosons was +// re-decayed hadronically. (Note that LHE format version larger than 3 may not +// require this hook). +// +// The hook inherits from PowhegHooks and as such it indirectly implements the +// following: +// - doVetoMPIStep +// - doVetoISREmission +// - doVetoMPIEmission +// and it overloads: +// - doVetoFSREmission, which works as follows (if POWHEG:veto==1): +// - if inResonance==true it vetoes all the emission that is harder than +// the scale of its parent (anti-)top quark or W^+(-) +// - if inResonance==false, it calls PowhegHooks::doVetoISREmission +// and it also implements: +// - doVetoProcessLevel, which is never used for vetoing (i.e. it always +// returns false). Instead it is used for the calculation of reconance scales +// using LHE kinematics. +// +// This version of the hooks is only suitable for use with fully compatible +// POWHEG-BOX Les Houches readers (such the one in main-PYTHIA82-lhef but +// not the one in main31.cc.) +// +// +// # Basic use +// +// In order to use this hook you must replace all the declarations and +// constructor calls of PowhegHooks to PowhegHooksBB4L: +// +// PowhegHooks *powhegHooks; -> PowhegHooksBB4L *powhegHooks; +// *powhegHooks = new PowhegHooks(); -> *powhegHooks = new PowhegHooksBB4L(); +// +// In order to switch it on set POWHEG:veto = 1 and +// POWHEG:bb4l:FSREmission:veto = 1. This will lead to a veto in ISR, FSR and +// MPI steps of pythia as expected using PowhegHooks in all the cases other than +// the case of FSR emission in resonance decay. Within resonance decays +// PowhegHooksBB4L takes over the control. +// +// Furthermore, this hook can also be used standalone without PowhegHooks, i.e. +// only the FSR emission from resonances will be vetoed (the default use in +// 1801.03944 and 1906.09166). In order to do that set +// POWHEG:bb4l:FSREmission:veto = 1 and POWHEG:veto = 0. Note that the this is +// not recommended for comparing against data but because it is easier to +// interpret it is often done in theoretical studies. +// +// Note that this version of the hook relies on the event "radtype" (1 for +// btilde, 2 for remnant) to be set by an external program, such as +// main-PYTHIA82-lhef in the radtype_ common block. +// There also exists a version of this hook in which the event "radtype" is +// read in from the .lhe file using pythia built in functionality. You need +// that version if you want to use this hook with main31.cc. +// +// +// # Expert use +// +// This hook also implements an alternative veto procedure which allows to +// assign a "SCALUP" type of scale to a resonance using the scaleResonance +// method. This is a much simpler veto but it is also clearly inferior as +// compared to the one implemented using the doVetoFSREmission method because +// the definition of the scale of the emission does not match the +// POWHEG-BOX-RES definition. Nevertheless, it can be activated using +// POWHEG:bb4l:ScaleResonance:veto = 1. Additionally one MUST switch off the +// other veto by calling on POWHEG:bb4l:FSREmission:veto = 0. +// +// The following settings are at the disposal of the user to control the +// behaviour of the hook +// - On/off switches for the veto: +// - POWHEG:bb4l:FSREmission:veto +// on/off switch for the default veto based on doFSREmission +// - POWHEG:bb4l:ScaleResonance:veto +// on/off switch for the alternative veto based on scaleResonance (only +// for expert users) +// - Important settings: +// - POWHEG:bb4l:ptMinVeto: MUST be set to the same value as the +// corresponding flag in POWHEG-BOX-RES +// - Alternatives for scale calculations +// - default: emission scale is calculated using POWHEG definitions and in +// the resonance rest frame +// - POWHEG:bb4l:FSREmission:vetoDipoleFrame: emission scale is calculated +// using POWHEG definitions in the dipole frame +// - POWHEG:bb4l:FSREmission:pTpythiaVeto: emission scale is calculated +// using Pythia definitions +// - Other flags: +// - POWHEG:bb4l:FSREmission:vetoQED: decides whether or not QED emission +// off quarks should also be vetoed (not implemented in the case of +// the ScaleResonance:veto) +// - POWHEG:bb4l:DEBUG: enables debug printouts on standard output #ifndef Pythia8_PowhegHooksBB4L_H #define Pythia8_PowhegHooksBB4L_H -// Includes #include "Pythia8/Pythia.h" #include -struct { - int radtype; -} radtype_; namespace Pythia8 { class PowhegHooksBB4L : public UserHooks { public: - //--- Constructor and destructor ------------------------------------------- - PowhegHooksBB4L() : nFSRvetoBB4l(0) {} - ~PowhegHooksBB4L() override { std::cout << "Number of FSR vetoed in BB4l = " << nFSRvetoBB4l << std::endl; } + PowhegHooksBB4L() {} + ~PowhegHooksBB4L() { std::cout << "Number of FSR vetoed in BB4l = " << nInResonanceFSRveto << std::endl; } - //--- Initialization ----------------------------------------------------------------------- - bool initAfterBeams() override { - // settings of this class + //--- Initialization ------------------------------------------------------- + bool initAfterBeams() { + // initialize settings of this class vetoFSREmission = settingsPtr->flag("POWHEG:bb4l:FSREmission:veto"); - onlyDistance1 = settingsPtr->flag("POWHEG:bb4l:FSREmission:onlyDistance1"); - dryRunFSR = settingsPtr->flag("POWHEG:bb4l:FSREmission:dryRun"); - vetoAtPL = settingsPtr->flag("POWHEG:bb4l:FSREmission:vetoAtPL"); - vetoQED = settingsPtr->flag("POWHEG:bb4l:FSREmission:vetoQED"); - vetoPartonLevel = settingsPtr->flag("POWHEG:bb4l:PartonLevel:veto"); - excludeFSRConflicting = settingsPtr->flag("POWHEG:bb4l:PartonLevel:excludeFSRConflicting"); - debug = settingsPtr->flag("POWHEG:bb4l:DEBUG"); - scaleResonanceVeto = settingsPtr->flag("POWHEG:bb4l:ScaleResonance:veto"); vetoDipoleFrame = settingsPtr->flag("POWHEG:bb4l:FSREmission:vetoDipoleFrame"); pTpythiaVeto = settingsPtr->flag("POWHEG:bb4l:FSREmission:pTpythiaVeto"); - //vetoProduction = (settingsPtr->mode("POWHEG:veto")==1); + vetoQED = settingsPtr->flag("POWHEG:bb4l:FSREmission:vetoQED"); + scaleResonanceVeto = settingsPtr->flag("POWHEG:bb4l:ScaleResonance:veto"); + debug = settingsPtr->flag("POWHEG:bb4l:DEBUG"); pTmin = settingsPtr->parm("POWHEG:bb4l:pTminVeto"); + vetoAllRadtypes = settingsPtr->flag("POWHEG:bb4l:vetoAllRadtypes"); + nInResonanceFSRveto = 0; return true; } //--- PROCESS LEVEL HOOK --------------------------------------------------- - - // called at the LHE level - inline bool canVetoProcessLevel() override { return true; } - inline bool doVetoProcessLevel(Event &e) override { + // This hook gets triggered for each event before the shower starts, i.e. at + // the LHE level. We use it to calculate the scales of resonances. + inline bool canVetoProcessLevel() { return true; } + inline bool doVetoProcessLevel(Event &e) { // extract the radtype from the event comment stringstream ss; // use eventattribute as comments not filled when using edm input - //ss << infoPtr->getEventComments(); ss << infoPtr->getEventAttribute("#rwgt"); string temp; - ss >> temp >> radtype_.radtype; + ss >> temp >> radtype; assert(temp == "#rwgt"); - + // we only calculate resonance scales for btilde events (radtype==1) + // remnant events are not vetoed + if (!vetoAllRadtypes && radtype == 2) + return false; // find last top and the last anti-top in the record - int i_top = -1, i_atop = -1; + int i_top = -1, i_atop = -1, i_wp = -1, i_wm = -1; for (int i = 0; i < e.size(); i++) { if (e[i].id() == 6) i_top = i; if (e[i].id() == -6) i_atop = i; + if (e[i].id() == 24) + i_wp = i; + if (e[i].id() == -24) + i_wm = i; } - if (i_top != -1) - topresscale = findresscale(i_top, e); - else - topresscale = 1e30; - if (i_top != -1) - atopresscale = findresscale(i_atop, e); - else - atopresscale = 1e30; - // initialize stuff - doVetoFSRInit(); - // do not veto, ever - return false; - } - - //--- PARTON LEVEL HOOK ---------------------------------------------------- + // if found calculate the resonance scale + topresscale = findresscale(i_top, e); + // similarly for anti-top + atopresscale = findresscale(i_atop, e); + // and for W^+ and W^- + wpresscale = findresscale(i_wp, e); + wmresscale = findresscale(i_wm, e); - // called after shower - bool retryPartonLevel() override { return vetoPartonLevel || vetoAtPL; } - inline bool canVetoPartonLevel() override { return vetoPartonLevel || vetoAtPL; } - inline bool doVetoPartonLevel(const Event &e) override { - if (radtype_.radtype == 2) - return false; - if (debug) { - if (dryRunFSR && wouldVetoFsr) { - double scale = getdechardness(vetoTopCharge, e); - cout << "FSRdecScale = " << vetoDecScale << ", PLdecScale = " << scale << ", ratio = " << vetoDecScale / scale - << endl; - } - } - if (vetoPartonLevel) { - double topdecscale = getdechardness(1, e); - double atopdecscale = getdechardness(-1, e); - if ((topdecscale > topresscale) || (atopdecscale > atopresscale)) { - //if(dryRunFSR && ! wouldVetoFsr) mydatacontainer_.excludeEvent = excludeFSRConflicting?1:0; - return true; - } else - //if(dryRunFSR && wouldVetoFsr) mydatacontainer_.excludeEvent = excludeFSRConflicting?1:0; - return false; - } - if (vetoAtPL) { - if (dryRunFSR && wouldVetoFsr) - return true; - else - return false; - } + // do not veto, ever return false; } //--- FSR EMISSION LEVEL HOOK ---------------------------------------------- - - // FSR veto: this should be true if we want PowhegHooksBB4l veto in decay - // OR PowhegHooks veto in production. (The virtual method - // PowhegHooks::canVetoFSREmission has been replaced by - // PowhegHooksBB4L::canVetoFSREmission, so FSR veto in production - // must be handled here. ISR and MPI veto are instead still - // handled by PowhegHooks.) - inline bool canVetoFSREmission() override { return vetoFSREmission; } // || vetoProduction; } - inline bool doVetoFSREmission(int sizeOld, const Event &e, int iSys, bool inResonance) override { - ////////////////////////////// - //VETO INSIDE THE RESONANCE // - ////////////////////////////// + // This hook gets triggered everytime the parton shower attempts to attach + // a FSR emission. + inline bool canVetoFSREmission() { return vetoFSREmission; } + inline bool doVetoFSREmission(int sizeOld, const Event &e, int iSys, bool inResonance) { + // FSR VETO INSIDE THE RESONANCE (if it is switched on) if (inResonance && vetoFSREmission) { + // get the participants of the splitting: the recoiler, the radiator and the emitted int iRecAft = e.size() - 1; int iEmt = e.size() - 2; int iRadAft = e.size() - 3; int iRadBef = e[iEmt].mother1(); - // find the top resonance the radiator originates from - int iTop = e[iRadBef].mother1(); - int distance = 1; - while (abs(e[iTop].id()) != 6 && iTop > 0) { - iTop = e[iTop].mother1(); - distance++; + // find the resonance the radiator originates from + int iRes = e[iRadBef].mother1(); + while (iRes > 0 && (abs(e[iRes].id()) != 6 && abs(e[iRes].id()) != 24)) { + iRes = e[iRes].mother1(); } - if (iTop == 0) { + if (iRes == 0) { infoPtr->errorMsg( - "Warning in PowhegHooksBB4L::doVetoFSREmission: emission in resonance not from top quark, not vetoing"); - return doVetoFSR(false, 0, 0); - //return false; + "Warning in PowhegHooksBB4L::doVetoFSREmission: emission in resonance not from the top quark or from the " + "W boson, not vetoing"); + return doVetoFSR(false, 0); } - int iTopCharge = (e[iTop].id() > 0) ? 1 : -1; + int iResId = e[iRes].id(); // calculate the scale of the emission double scale; @@ -151,132 +198,121 @@ namespace Pythia8 { scale = pTpythia(e, iRadAft, iEmt, iRecAft); //.. or using POWHEG pT definition else { - Vec4 pr(e[iRadAft].p()), pe(e[iEmt].p()), pt(e[iTop].p()), prec(e[iRecAft].p()), psystem; + Vec4 pr(e[iRadAft].p()), pe(e[iEmt].p()), pres(e[iRes].p()), prec(e[iRecAft].p()), psystem; // The computation of the POWHEG pT can be done in the top rest frame or in the diple one. // pdipole = pemt +prec +prad (after the emission) // For the first emission off the top resonance pdipole = pw +pb (before the emission) = ptop if (vetoDipoleFrame) psystem = pr + pe + prec; else - psystem = pt; + psystem = pres; // gluon splitting into two partons if (e[iRadBef].id() == 21) scale = gSplittingScale(psystem, pr, pe); // quark emitting a gluon (or a photon) - else if (abs(e[iRadBef].id()) == 5 && ((e[iEmt].id() == 21) && !vetoQED)) + else if (abs(e[iRadBef].id()) <= 5 && ((e[iEmt].id() == 21) && !vetoQED)) scale = qSplittingScale(psystem, pr, pe); // other stuff (which we should not veto) - else + else { scale = 0; + } } - if (iTopCharge > 0) { - if (onlyDistance1) { - if (debug && (distance == 1) && scale > topresscale && !wouldVetoFsr) - cout << e[iTop].id() << ": " << e[iRadBef].id() << " > " << e[iRadAft].id() << " + " << e[iEmt].id() - << "; " << scale << endl; - return doVetoFSR((distance == 1) && scale > topresscale, scale, iTopCharge); - } else { - if (debug && scale > topresscale && !wouldVetoFsr) - cout << e[iTop].id() << ": " << e[iRadBef].id() << " > " << e[iRadAft].id() << " + " << e[iEmt].id() - << "; " << scale << endl; - return doVetoFSR(scale > topresscale, scale, iTopCharge); - } - } else if (iTopCharge < 0) { - if (onlyDistance1) { - if (debug && (distance == 1) && scale > atopresscale && !wouldVetoFsr) - cout << e[iTop].id() << ": " << e[iRadBef].id() << " > " << e[iRadAft].id() << " + " << e[iEmt].id() - << "; " << scale << endl; - return doVetoFSR((distance == 1) && scale > atopresscale, scale, iTopCharge); - } else { - if (debug && scale > topresscale && !wouldVetoFsr) - cout << e[iTop].id() << ": " << e[iRadBef].id() << " > " << e[iRadAft].id() << " + " << e[iEmt].id() - << "; " << scale << endl; - return doVetoFSR(scale > atopresscale, scale, iTopCharge); - } + // compare the current splitting scale to the correct resonance scale + if (iResId == 6) { + if (debug && scale > topresscale) + cout << iResId << ": " << e[iRadBef].id() << " > " << e[iRadAft].id() << " + " << e[iEmt].id() << "; " + << scale << endl; + return doVetoFSR(scale > topresscale, scale); + } else if (iResId == -6) { + if (debug && scale > atopresscale) + cout << iResId << ": " << e[iRadBef].id() << " > " << e[iRadAft].id() << " + " << e[iEmt].id() << "; " + << scale << endl; + return doVetoFSR(scale > atopresscale, scale); + } else if (iResId == 24) { + if (debug && scale > wpresscale) + cout << iResId << ": " << e[iRadBef].id() << " > " << e[iRadAft].id() << " + " << e[iEmt].id() << "; " + << scale << endl; + return doVetoFSR(scale > wpresscale, scale); + } else if (iResId == -24) { + if (debug && scale > wmresscale) + cout << iResId << ": " << e[iRadBef].id() << " > " << e[iRadAft].id() << " + " << e[iEmt].id() << "; " + << scale << endl; + return doVetoFSR(scale > wmresscale, scale); } else { - cout << "Bug in PohwgeHooksBB4l" << endl; + infoPtr->errorMsg("Error in PowhegHooksBB4L::doVetoFSREmission: unimplemented case"); + exit(-1); } } - ///////////////////////////////// - // VETO THE PRODUCTION PROCESS // - ///////////////////////////////// - // covered by multiuserhook, i.e. need to turn on EV1 - // else if(!inResonance && vetoProduction){ - // return EmissionVetoHook1::doVetoFSREmission(sizeOld, e, iSys, inResonance); - // } - - return false; + // In CMSSW, the production process veto is done in EmissionVetoHook1.cc + // so for events outside resonance, nothing needs to be done here + else { + return false; + } } - inline bool doVetoFSR(bool condition, double scale, int iTopCharge) { - if (radtype_.radtype == 2) + inline bool doVetoFSR(bool condition, double scale) { + if (!vetoAllRadtypes && radtype == 2) return false; if (condition) { - if (!wouldVetoFsr) { - wouldVetoFsr = true; - vetoDecScale = scale; - vetoTopCharge = iTopCharge; - } - if (dryRunFSR) - return false; - else { - nFSRvetoBB4l++; - return true; - } - } else - return false; - } - - inline void doVetoFSRInit() { - wouldVetoFsr = false; - vetoDecScale = -1; - vetoTopCharge = 0; + nInResonanceFSRveto++; + return true; + } + return false; } //--- SCALE RESONANCE HOOK ------------------------------------------------- // called before each resonance decay shower - inline bool canSetResonanceScale() override { return scaleResonanceVeto; } - // if the resonance is the (anti)top set the scale to: - // ---> (anti)top virtuality if radtype=2 - // ---> (a)topresscale otherwise + inline bool canSetResonanceScale() { return scaleResonanceVeto; } + // if the resonance is the (anti)top or W+/W- set the scale to: + // - if radtype=2 (remnant): resonance virtuality + // - if radtype=1 (btilde): + // - (a)topresscale/wp(m)resscale for tops and Ws + // - a large number otherwise // if is not the top, set it to a big number - inline double scaleResonance(int iRes, const Event &e) override { - if (e[iRes].id() == 6) { - if (radtype_.radtype == 2) - return sqrt(e[iRes].m2Calc()); - else + inline double scaleResonance(int iRes, const Event &e) { + if (!vetoAllRadtypes && radtype == 2) + return sqrt(e[iRes].m2Calc()); + else { + if (e[iRes].id() == 6) return topresscale; - } else if (e[iRes].id() == -6) { - if (radtype_.radtype == 2) - return sqrt(e[iRes].m2Calc()); - else + else if (e[iRes].id() == -6) return atopresscale; - } else - return pow(10.0, 30.); + else if (e[iRes].id() == 24) + return wpresscale; + else if (e[iRes].id() == 24) + return wmresscale; + else + return 1e30; + } } //--- Internal helper functions -------------------------------------------- - // Calculates the scale of the hardest emission from within the resonance system // translated by Markus Seidel modified by Tomas Jezo inline double findresscale(const int iRes, const Event &event) { - double scale = 0.; + // return large scale if the resonance position is ill defined + if (iRes < 0) + return 1e30; + // get number of resonance decay products int nDau = event[iRes].daughterList().size(); + // iRes is not decayed, return high scale equivalent to + // unrestricted shower if (nDau == 0) { - // No resonance found, set scale to high value - // Pythia will shower any MC generated resonance unrestricted - scale = 1e30; - } else if (nDau < 3) { - // No radiating resonance found - scale = pTmin; - } else if (abs(event[iRes].id()) == 6) { - // Find top daughters + return 1e30; + } + // iRes did not radiate, this means that POWHEG pt scale has + // evolved all the way down to pTmin + else if (nDau < 3) { + return pTmin; + } + // iRes is a (anti-)top quark + else if (abs(event[iRes].id()) == 6) { + // find top daughters int idw = -1, idb = -1, idg = -1; - for (int i = 0; i < nDau; i++) { int iDau = event[iRes].daughterList()[i]; if (abs(event[iDau].id()) == 24) @@ -290,20 +326,48 @@ namespace Pythia8 { // Get daughter 4-vectors in resonance frame Vec4 pw(event[idw].p()); pw.bstback(event[iRes].p()); - Vec4 pb(event[idb].p()); pb.bstback(event[iRes].p()); + Vec4 pg(event[idg].p()); + pg.bstback(event[iRes].p()); + + // Calculate scale and return it + return sqrt(2 * pg * pb * pg.e() / pb.e()); + } + // iRes is a W+(-) boson + else if (abs(event[iRes].id()) == 24) { + // Find W daughters + int idq = -1, ida = -1, idg = -1; + for (int i = 0; i < nDau; i++) { + int iDau = event[iRes].daughterList()[i]; + if (event[iDau].id() == 21) + idg = iDau; + else if (event[iDau].id() > 0) + idq = iDau; + else if (event[iDau].id() < 0) + ida = iDau; + } + // Get daughter 4-vectors in resonance frame + Vec4 pq(event[idq].p()); + pq.bstback(event[iRes].p()); + Vec4 pa(event[ida].p()); + pa.bstback(event[iRes].p()); Vec4 pg(event[idg].p()); pg.bstback(event[iRes].p()); // Calculate scale - scale = sqrt(2 * pg * pb * pg.e() / pb.e()); - } else { - scale = 1e30; + Vec4 pw = pq + pa + pg; + double q2 = pw * pw; + double csi = 2 * pg.e() / sqrt(q2); + double yq = 1 - pg * pq / (pg.e() * pq.e()); + double ya = 1 - pg * pa / (pg.e() * pa.e()); + // and return it + return sqrt(min(1 - yq, 1 - ya) * pow2(csi) * q2 / 2); } - - return scale; + // in any other case just return a high scale equivalent to + // unrestricted shower + return 1e30; } // The following routine will match daughters of particle `e[iparticle]` to an expected pattern specified via the list of expected particle PDG ID's `ids`, @@ -393,109 +457,27 @@ namespace Pythia8 { return (sqrt(pTnow)); } - inline double getdechardness(int topcharge, const Event &e) { - int tid = 6 * topcharge, wid = 24 * topcharge, bid = 5 * topcharge, gid = 21, wildcard = 0; - // find last top in the record - int i_top = -1; - Vec4 p_top, p_b, p_g, p_g1, p_g2; - for (int i = 0; i < e.size(); i++) - if (e[i].id() == tid) { - i_top = i; - p_top = e[i].p(); - } - if (i_top == -1) - return -1.0; - - // summary of cases - // 1.) t > W b - // a.) b > 3 ... error - // b.) b > b g ... h = sqrt(2*p_g*p_b*p_g.e()/p_b.e()) - // c.) b > other ... h = -1 - // return h - // 2.) t > W b g - // a.) b > 3 ... error - // b.) b > b g ... h1 = sqrt(2*p_g*p_b*p_g.e()/p_b.e()) - // c.) b > other ... h1 = -1 - // i.) g > 3 ... error - // ii.) g > 2 ... h2 = sqrt(2*p_g1*p_g2*p_g1.e()*p_g2.e()/(pow(p_g1.e(),2)+pow(p_g2.e(),2))) ); - // iii.) g > other ... h2 = -1 - // return max(h1,h2) - // 3.) else ... error - - vector momenta; - vector positions; - - // 1.) t > b W - if (match_decay(i_top, e, vector{wid, bid}, positions, momenta, false)) { - double h; - int i_b = positions[1]; - // a.+b.) b > 3 or b > b g - if (match_decay(i_b, e, vector{bid, gid}, positions, momenta)) - h = qSplittingScale(e[i_top].p(), momenta[0], momenta[1]); - // c.) b > other - else - h = -1; - return h; - } - // 2.) t > b W g - else if (match_decay(i_top, e, vector{wid, bid, gid}, positions, momenta, false)) { - double h1, h2; - int i_b = positions[1], i_g = positions[2]; - // a.+b.) b > 3 or b > b g - if (match_decay(i_b, e, vector{bid, gid}, positions, momenta)) - h1 = qSplittingScale(e[i_top].p(), momenta[0], momenta[1]); - // c.) b > other - else - h1 = -1; - // i.+ii.) g > 3 or g > 2 - if (match_decay(i_g, e, vector{wildcard, wildcard}, positions, momenta)) - h2 = gSplittingScale(e[i_top].p(), momenta[0], momenta[1]); - // c.) b > other - else - h2 = -1; - return max(h1, h2); - } - // 3.) else - else { - cout << "getdechardness" << endl; - cout << "top at position " << i_top << endl; - cout << "with " << e[i_top].daughterList().size() << " daughters " << endl; - for (unsigned i = 0; i < e[i_top].daughterList().size(); i++) { - int di = e[i_top].daughterList()[i]; - cout << "with daughter " << di << ": " << e[di].id() << endl; - } - exit(-1); - } - } - - //-------------------------------------------------------------------------- - - // Functions to return information - - // inline int getNFSRveto() { return nFSRveto; } + // Functions to return statistics about the veto + inline int getNInResonanceFSRVeto() { return nInResonanceFSRveto; } //-------------------------------------------------------------------------- private: // FSR emission veto flags - bool vetoFSREmission, dryRunFSR, wouldVetoFsr, onlyDistance1, vetoAtPL, vetoQED; - // Parton Level veto flags - bool vetoPartonLevel, excludeFSRConflicting; - // Scale Resonance veto flags + bool vetoFSREmission, vetoQED; + // scale Resonance veto flags double scaleResonanceVeto; // other flags bool debug; - // internal: resonance scales - double topresscale, atopresscale; - // internal: inter veto communication - double vetoDecScale; - int vetoTopCharge; bool vetoDipoleFrame; bool pTpythiaVeto; - //bool vetoProduction; double pTmin; - // Statistics on vetos - unsigned long int nFSRvetoBB4l; + bool vetoAllRadtypes; + // veto counter + int nInResonanceFSRveto; + // internal: resonance scales + double topresscale, atopresscale, wpresscale, wmresscale; + int radtype; }; //========================================================================== diff --git a/GeneratorInterface/Pythia8Interface/src/Py8InterfaceBase.cc b/GeneratorInterface/Pythia8Interface/src/Py8InterfaceBase.cc index 99f8ec98dfcd1..885a69a427c33 100644 --- a/GeneratorInterface/Pythia8Interface/src/Py8InterfaceBase.cc +++ b/GeneratorInterface/Pythia8Interface/src/Py8InterfaceBase.cc @@ -112,18 +112,14 @@ namespace gen { //add settings for powheg resonance scale calculation fMasterGen->settings.addFlag("POWHEGres:calcScales", false); fMasterGen->settings.addFlag("POWHEG:bb4l", false); - fMasterGen->settings.addFlag("POWHEG:bb4l:FSREmission:onlyDistance1", false); fMasterGen->settings.addFlag("POWHEG:bb4l:FSREmission:veto", false); - fMasterGen->settings.addFlag("POWHEG:bb4l:FSREmission:dryRun", false); - fMasterGen->settings.addFlag("POWHEG:bb4l:FSREmission:vetoAtPL", false); fMasterGen->settings.addFlag("POWHEG:bb4l:FSREmission:vetoQED", false); - fMasterGen->settings.addFlag("POWHEG:bb4l:PartonLevel:veto", false); - fMasterGen->settings.addFlag("POWHEG:bb4l:PartonLevel:excludeFSRConflicting", false); fMasterGen->settings.addFlag("POWHEG:bb4l:DEBUG", false); fMasterGen->settings.addFlag("POWHEG:bb4l:ScaleResonance:veto", false); fMasterGen->settings.addFlag("POWHEG:bb4l:FSREmission:vetoDipoleFrame", false); fMasterGen->settings.addFlag("POWHEG:bb4l:FSREmission:pTpythiaVeto", false); fMasterGen->settings.addParm("POWHEG:bb4l:pTminVeto", 10.0, true, true, 0.0, 10.); + fMasterGen->settings.addFlag("POWHEG:bb4l:vetoAllRadtypes", false); fMasterGen->setRndmEnginePtr(&p8RndmEngine_); fDecayer->setRndmEnginePtr(&p8RndmEngine_);