Drastic parton showering effects at NLO
Hi there,
I am currently running VBF Higgs production in a BSM model at NLO with Pythia8. However, after parton showering, the \Delta \phi of leading and subleading jet distribution has been changed drastically. After showering, it is dominated at \Delta \phi = \pi, which is always backtoback in xy plane, however, it is not necessary to be like this, and it is quite different from the \Delta \phi distribution in partionlevel (no PS).
Is that possible that some potential problems occur when passing NLO computation to Pythia8 inside MadGraph5 interface? Or somewhere I need to be careful when I use MadGraph and Pythia8 at NLO?
Thank you very much!
Question information
 Language:
 English Edit question
 Status:
 Solved
 Assignee:
 Paolo Torrielli Edit question
 Solved by:
 Paolo Torrielli
 Solved:
 Last query:
 Last reply:
Revision history for this message

#1 
Dear Hantian,
for Dphi I'd normally expect way more events at Dphi=pi than at Dphi=0,
both for fixed NLO and for NLO+PS, with shower radiation tending to
flatten the distribution (i.e. to fill the Dphi=0 region) with respect to fixed
order. Without looking at plots it is a bit difficult to have an idea of the
drastic difference you see, though.
If you are asking about potential issues linked to the use of the default
Hjj analysis provided with the package, we have used it many times (in
the SM VBF) without noticing suspect features. If you are asking something
else, could you kindly specify what you mean by 'potential problems'?
In any case, a LO+PS run with the same setup, and/or a (N)LO+PS run with
a different shower (e.g. HW++) would be beneficial in order to understand
the origin of these features.
Cheers.
Paolo
Revision history for this message

#2 
Hi Paolo,
Thank you for your message. However, I checked the VBF Higgs production in SM model by default, the showering effects of Pythia8 does not flatten the distribution. More precisely, the no PS distribution is rather flat, with a bit more events on \Delta \phi = \pi than \Delta \phi =0 (this confirms the results of Maltoni's paper ArXiv: 1311.1829), but after showering the distribution is drastically changed to have pick at \Delta \phi = \pi with very very few events on \Delta \phi = 0.
Then I would guess it is maybe the problem associated to parton shower program...
Best,
Hantian
Revision history for this message

#3 
Hi Hantian,
OK: from the paper you cited, fig. 2, the effect of parton showering over
fixed NLO is rather minimal for Dphi (even smaller than I remembered,
thanks for pointing out that reference).
Are you sure you’re not clustering the decay products of the resonance
into jets at NLO+PS? If this were the case, in the fixed NLO and in the
NLO+PS simulations you would be considering very different jets, that
could explain a dramatic effect.
Moreover, what happens at LO+PS? Do you also see a sizeable distortion?
Thanks.
Cheers.
Paolo
> On 14 Aug 2017, at 17:57, Hantian Zhang <email address hidden> wrote:
>
> Question #655841 on MadGraph5_aMC@NLO changed:
> https:/
>
> Hantian Zhang posted a new comment:
> Hi Paolo,
>
> Thank you for your message. However, I checked the VBF Higgs production
> in SM model by default, the showering effects of Pythia8 does not
> flatten the distribution. More precisely, the no PS distribution is
> rather flat, with a bit more events on \Delta \phi = \pi than \Delta
> \phi =0 (this confirms the results of Maltoni's paper ArXiv: 1311.1829),
> but after showering the distribution is drastically changed to have pick
> at \Delta \phi = \pi with very very few events on \Delta \phi = 0.
>
> Then I would guess it is maybe the problem associated to parton shower
> program...
>
> Best,
> Hantian
>
> 
> You received this question notification because you are assigned to this
> question.
Revision history for this message

#4 
Hi Paolo,
I checked both NLO+PS and LO+PS, PS all have the same drastic effects on Dphi distribution.
For the clustering, I use FastJet in MadAnalysis5 with my own Dphi observable (since MA5 Dphi is not correct in current version).
Would you like to provide more information regards the "clustering the decay products of the resonance into jets at NLO+PS"? Then I could look into details on my simulation to have a check.
Thank you very much!
Best,
Hantian
Revision history for this message

#5 
Hi Hantian,
> Would you like to provide more information regards the "clustering the
> decay products of the resonance into jets at NLO+PS"? Then I could look
> into details on my simulation to have a check.
What I mean is: if, in the showered simulation, the decay chain of
your (Higgslike) resonance ends up with jets, then the shape of all
jet observables will be very much distorted with respect to fixed order,
where (I guess) you keep the resonance stable.
But this distortion would come only from having (artificially) no jets
from the decay in the fixedorder calculation.
Are you keeping the Higgs stable in your PS simulation (or are you
decaying it to objects that are not clustered into jets)?
Thanks.
Cheers.
Paolo
Revision history for this message

#6 
Hi Paolo,
Thank you for pointing this. I didn't keep Higgs stable in PS before, which is set to h_stable = F by default. However, after I setting it to be h_stable = T, the drastic PS effects are still there.
Best,
Hantian
Revision history for this message

#7 
Hi Hantian,
OK, strange, I would have thought of an effect of this kind.
Just to be sure, please double check that the events after
showering indeed have stable Higgses (just to make sure that
command in the shower card is correctly understood. Alternatively
you could decay H > tau tau in the card.).
I’ll think about other sources of this effect. Could you maybe
paste your shower card and run card?
Thanks.
Cheers.
Paolo
> On 14 Aug 2017, at 20:03, Hantian Zhang <email address hidden> wrote:
>
> Question #655841 on MadGraph5_aMC@NLO changed:
> https:/
>
> Hantian Zhang posted a new comment:
> Hi Paolo,
>
> Thank you for pointing this. I didn't keep Higgs stable in PS before,
> which is set to h_stable = F by default. However, after I setting it to
> be h_stable = T, the drastic PS effects are still there.
>
> Best,
> Hantian
>
> 
> You received this question notification because you are assigned to this
> question.
Revision history for this message

#8 
Hi Paolo,
Thank you very much! Here is the shower card and run card:
#******
# MadGraph5_aMC@NLO *
# *
# shower_card.dat aMC@NLO *
# *
# This file is used to set the parameters for the shower. *
# *
# Some notation/
# *
# Lines starting with a hash (#) are info or comments *
# *
# mind the format: variable = value # comment *
#******
#
#******
# Shower settings *
#******
# Number of events, jobs, errors, and random seeds *
#******
nevents = 1 # N evts to shower (< 0 = all)
nsplit_jobs = 1 # N jobs to run in parallel (< 100!!)
combine_td = T # combine the topdrawer/HwU files if nsplit_jobs>1
maxprint = 2 # N evts to print in the log
maxerrs = 0.1 # max fraction of errors
rnd_seed = 0 # 1st random seed (0 = default)
rnd_seed2 = 0 # 2nd random seed (0 = default) !ONLY FOR HWERIG6!
#******
# PDFs and nonperturbative modelling *
#******
pdfcode = 0 # 0 = internal, 1 = same as NLO, other = lhaglue
ue_enabled = F # underlying event
hadronize = T # hadronisation on/off !IGNORED BY HERWIG6!
lambda_5 = 1 # Lambda_5 (< 0 = default) !IGNORED BY PYTHIA8!
#******
# Stable or unstable particles *
#******
b_stable = F # set B hadrons stable
pi_stable = T # set pi0's stable
wp_stable = F # set w+'s stable
wm_stable = F # set w's stable
z_stable = F # set z0's stable
h_stable = T # set Higgs' stable
tap_stable = F # set tau+'s stable
tam_stable = F # set tau's stable
mup_stable = F # set mu+'s stable
mum_stable = F # set mu's stable
#******
# Mass of the b quark *
#******
b_mass = 1 # if < 0 = read from SubProcesses/
#******
# Special settings *
#******
is_4lep = F # T if 4lepton production !ONLY FOR PYTHIA6!
is_bbar = F # T if bb~ production !ONLY FOR HERWIG6!
#******
# FxFx merging parameters !ONLY FOR PYTHIA8!
#******
Qcut = 1.0 # Merging scale
njmax = 0 # Maximal multiplicity in the merging
#******
# Decay channels *
#******
# Syntax for HERWIG6 *
# DM_I = M > D1 D2 @ BR @ ME *
# corresponding to call to HWMODK(
# I < 100, M is the decaying resonance, D1, D2, ... are the decay *
# products (up to five), BR is the branching ratio and ME is the type *
# of matrix element to be used in the decay. *
# BR's are correctly understood only if they add up to 1, and only if *
# no more than three modes are required for a given resonance. *
# WARNING: the order of decay products in > 2body decays IS RELEVANT. *
# *
# Syntax for PYTHIA6 *
# DM_I = M > D1 D2 @ BR @ ME *
# WARNING: turning hadronisation off disables top decays *
# WARNING: 1 > n decays (with n > 2) are handled through a sequence *
# of 1 > 2 decays. *
# WARNING: entries BR and ME are ignored *
# *
# Syntax for HERWIG++ *
# DM_I = M > D1 D2 @ BR @ ME *
# WARNING: entries BR and ME are ignored *
# *
# Syntax for PYTHIA8 *
# DM_I = M:onIfAny = D1 D2 *
# or similar, according to the offical PYTHIA8 decay syntax, see *
# the online PYTHIA8 manual *
# WARNING: 1 > n decays (with n > 2) are handled through a sequence *
# of 1 > 2 decays. *
# *
# Examples *
# Z > e+ e or mu+ mu with BR = 0.5 each, HERWIG6 *
# DM_1 = 23 > 11 11 @ 0.5d0 @ 100
# DM_2 = 23 > 13 13 @ 0.5d0 @ 100
# H > ta+ ta with BR = 1, HERWIG6 or HERWIG++ *
# DM_3 = 25 > 15 15 @ 1.0d0 @ 0
# t > ve e+ b with BR = 1, HERWIG6 or HERWIG++ *
# DM_4 = 6 > 12 11 5 @ 1d0 @ 100
# t > ve e+ b with BR = 1, PYTHIA6 *
# DM_5 = 6 > 24 5 @ 1d0 @ 100
# DM_6 = 24 > 12 11 @ 1d0 @ 100
# W+ > ve e+, W > vm~ mu, PYTHIA8 *
# DM_1 = 24:onMode = off
# DM_2 = 24:onPosIfAny = 11 12
# DM_3 = 24:onNegIfAny = 13 14
# W+ > ve e+, W > ve~ e and vm~ mu, PYTHIA8 *
# DM_1 = 24:onMode = off
# DM_2 = 24:onIfAny = 11 12
# DM_3 = 24:onNegIfAny = 13 14
#******
# Extra libraries/analyses *
#******
# The following lines need to be changed if the user does not want to *
# create a StdHEP/HepMC file, but to directly run an own analysis (to *
# be placed in HWAnalyzer or analogous MCatNLO subfolders). *
# Please use files in those folders as examples. *
# Remember that if your analysis uses hbook or is in the HwU format, *
# you must also add to hbook.o or HwU.o to the ANALYSE list as well. *
#******
extralibs = pythia8 # Extralibraries (not LHAPDF)
extrapaths = ../lib /app/cern/
INCLUDEPATHS = # Path to header files needed by c++
ANALYSE = # User's analysis and histogramming
#******
# MadGraph5_aMC@NLO *
# *
# run_card.dat aMC@NLO *
# *
# This file is used to set the parameters of the run. *
# *
# Some notation/
# *
# Lines starting with a hash (#) are info or comments *
# *
# mind the format: value = variable ! comment *
# *
# Some of the values of variables can be list. These can either be *
# comma or space separated. *
#******
#
#******
# Running parameters
#******
#
#******
# Tag name for the run (one word) *
#******
tag_1 = run_tag ! name of the run
#******
# Number of LHE events (and their normalization) and the required *
# (relative) accuracy on the Xsec. *
# These values are ignored for fixed order runs *
#******
50000 = nevents ! Number of unweighted events requested
1.0 = req_acc ! Required accuracy (1=auto determined from nevents)
1 = nevt_job! Max number of events per job in event generation.
! (1= no split).
#******
# Normalize the weights of LHE events such that they sum or average to *
# the total cross section *
#******
average = event_norm ! average or sum
#******
# Number of points per itegration channel (ignored for aMC@NLO runs) *
#******
0.01 = req_acc_FO ! Required accuracy (1=ignored, and use the
# These numbers are ignored except if req_acc_FO is equal to 1
5000 = npoints_FO_grid ! number of points to setup grids
4 = niters_FO_grid ! number of iter. to setup grids
10000 = npoints_FO ! number of points to compute Xsec
6 = niters_FO ! number of iter. to compute Xsec
#******
# Random number seed *
#******
0 = iseed ! rnd seed (0=assigned automatically=
#******
# Collider type and energy *
#******
1 = lpp1 ! beam 1 type (0 = no PDF)
1 = lpp2 ! beam 2 type (0 = no PDF)
6500.0 = ebeam1 ! beam 1 energy in GeV
6500.0 = ebeam2 ! beam 2 energy in GeV
#******
# PDF choice: this automatically fixes also alpha_s(MZ) and its evol. *
#******
nn23nlo = pdlabel ! PDF set
244600 = lhaid ! If pdlabel=lhapdf, this is the lhapdf number. Only
! numbers for central PDF sets are allowed. Can be a list;
! PDF sets beyond the first are included via reweighting.
#******
# Include the NLO Monte Carlo subtr. terms for the following parton *
# shower (HERWIG6  HERWIGPP  PYTHIA6Q  PYTHIA6PT  PYTHIA8) *
# WARNING: PYTHIA6PT works only for processes without FSR!!!! *
#******
PYTHIA8 = parton_shower
1.0 = shower_scale_factor ! multiply default shower starting
#******
# Renormalization and factorization scales *
# (Default functional form for the nonfixed scales is the sum of *
# the transverse masses divided by two of all final state particles *
# and partons. This can be changed in SubProcesses/
# dynamical_
#******
False = fixed_ren_scale ! if .true. use fixed ren scale
False = fixed_fac_scale ! if .true. use fixed fac scale
91.118 = muR_ref_fixed ! fixed ren reference scale
91.118 = muF_ref_fixed ! fixed fact reference scale
1 = dynamical_
! dynamical choices. Can be a list; scale choices beyond the
! first are included via reweighting
1.0 = muR_over_ref ! ratio of current muR over reference muR
1.0 = muF_over_ref ! ratio of current muF over reference muF
#******
# Reweight variables for scale dependence and PDF uncertainty *
#******
1.0, 2.0, 0.5 = rw_rscale ! muR factors to be included by reweighting
1.0, 2.0, 0.5 = rw_fscale ! muF factors to be included by reweighting
True = reweight_scale ! Reweight to get scale variation using the
! rw_rscale and rw_fscale factors. Should be a list of
! booleans of equal length to dynamical_
! specify for which choice to include scale dependence.
False = reweight_PDF ! Reweight to get PDF uncertainty. Should be a
! list booleans of equal length to lhaid to specify for
! which PDF set to include the uncertainties.
#******