What is SysCalc

=Running the code

Standalone Interface

scalefact	various value used for the variation of the factorization and renormalization scale
scalecorrelation	Define how the correlation between the factorization and renormalization scale has to be handle	-1: all combination
		-2 only correlated ones
alpsfact |
matchscale |
PDF |

This is an example of configuration file { # Central scale factors scalefact: 0.5 1 2 # Scale correlation # Special value -1: all combination (N2) # Special value -2: only correlated variation # Otherwise list of index N*fac_index + ren_index #index starts at 0 scalecorrelation: 0 # \alpha_s emission scale factors alpsfact: 0.5 1 2 # matching scales matchscale: 32 50 80 # PDF sets and number of members (optional) PDF: CT10.LHgrid 53 MSTW2008nlo68cl.LHgrid }

MadEvent Interface

In Order to use SysCalc from the MadEvent interface you need to

Technical Details

=======================================================
Description of variables stored for systematics studies
=======================================================
By: Johan Alwall, 28/3/2012

-------------------------------------------------------------------
Turn on systematics info with the flag use_syst in the run_card.dat
Note that systematics only works with matching (ickkw set to 1)
-------------------------------------------------------------------

Parameters that can be varied after-the-fact (without need to rerun
Pythia+detector simulation):

- Central renormalization scale
- Central factorization scale
- PDF choice
- Emission renormalization scale factor
- PDF reweighting scale factor (not available at present)
- QCUT scale

=============================================
Variational parameters in reweight.f:
=============================================

*********************************************
Central scale (ren scale):
*********************************************
line 702: scale (asref = alphas(scale))

Event weight given by:
----------------------
alpha_s weight: alphas(scale)^N where 
                N=#(QCD vertices) - #(emission alpha_s vertices below)

*********************************************
Emission alpha_s reweighting:
*********************************************
line 912: qnow = sqrt(q2now)

Event weight given by:
----------------------
alphas(alpsfact*qnow)

*********************************************
PDF reweighting:
*********************************************
line 873  (initial pdf): ipdgcl(idacl(n,i)), xnow(j), q2now (pdgini,xini, q2ini)
line 1070 (cont.   pdf): ipdgcl(idacl(n,i)), xnow(j), q2now (pdgint,xint, q2int)
                         etc. (for both sides (1,2))

Event weight given by:
----------------------
initial pdf(pdgini,xini,q2ini)
*pdf(pdgint,xint,q2int)/pdf(pdgint,xint,q2ini)
... etc.
Note: Central fact scale variation corresponds to reweighting the last 
      scale only on each side

=============================================
Variation of QCUT in ME2pythia.f:
=============================================

failing criteria in parentheses. line numbers approximate.

SHOWERKT: QCUT
line 999: PTSORT(1) (in lhe file) (< QCUT)
line 1012 (non-highest mult): shower kt (> QCUT)
line 1025 (highest mult): shower kt (> PTSORT(1))

kT-MLM: YCUT=QCUT**2
line 1090: NJETS (< NLJETS)
           actually Y(NLJETS) (< YCUT)
line 1107 (non-highest mult): Y(NLJETS+1) (> YCUT)
if highest mult case: YCUT=PTSORT(1)
line 1133: Y(NN) (> YCUT)
line 1145: If not clustered, fail
line 1176: Y(2) (> YCUT)

So, just need three number for systematic variation of QCUT:
SMIN. Fail if < QCUT: 
   For SHOWERKT: PTSORT(1)
   For kT-MLM: Y(NLJETS)
SCOMP. Comparison number: 
   For highest mult: max(QCUT,PTSORT(1))
   Otherwise QCUT
   Perhaps use minimum safe QCUT or 0 instead of QCUT for systematics studies
   In any case, use max(QCUTcurr,comparison number) for arbitrary QCUTcurr.
SMAX. Fail if > comparison number:
   For SHOWERKT: shower kt
   For kt-MLM: max(Y(NLJETS+1),Y(NN),Y(2))

Note that some events will always fail - I suggest to simply ignore
those (as well as requiring minimum safe QCUT = xqcut for SHOWERKT and
max(xqcut+10,xqcut*1.3) for kT-MLM).

Event weight given by:
----------------------
1 if QCUT < SMIN and SMAX < max(QCUT, SCOMP)
otherwise 0


===============================================
Each line in the syst.dat file has the entries:
===============================================

<mgrwt event="event_num">
<rscale>n_qcd ren_scale</rscale>
<asrwt>n_alpsem alpsem_scale(1) ... alpsem_scale(n_alpsem)</asrwt>
<pdfrwt beam="1">n_pdfrw1 pdf_pdg_code1(1) ... pdf_pdg_code1(n_pdgrw1) \
 pdf_x1(1) ... pdf_x1(n_pdfrw1) pdf_q1(1) ... pdf_q1(n_pdfrw1)</pdfwgt>
<pdfrwt beam="2">n_pdfrw2 pdf_pdg_code2(1) ... pdf_pdg_code2(n_pdgrw2) \
 pdf_x2(1) ... pdf_x2(n_pdfrw2) pdf_q2(1) ... pdf_q2(n_pdfrw2)</pdfwgt>
<totfact>total_reweight_factor</totfact>
<matchscale>SMIN SCOMP SMAX</matchscale>
</mgrwt>

Total event weight for event event_num given by:
------------------------------------------------
alpha_s(scalefact*ren_scale)^(n_qcd) *  # central ren scale
alpha_s(alpsfact*alpsem_scale(1)) *              # emission ren scale
alpha_s(alpsfact*alpsem_scale(2)) *              # emission ren scale
...
pdf(pdf_pdg_code1(1),pdf_x1(1),pdf_q1(1))*       # initial state pdf
pdf(pdf_pdg_code1(2),pdf_x1(2),pdf_q1(2))/
pdf(pdf_pdg_code1(2),pdf_x1(2),pdf_q1(2))*       # pdf reweighting
pdf(pdf_pdg_code1(3),pdf_x1(3),pdf_q1(3))/
pdf(pdf_pdg_code1(3),pdf_x1(3),pdf_q1(2))*       # pdf reweighting
 ...
pdf(pdf_pdg_code1(n_pdfrw1),pdf_x1(n_pdfrw1),scalefact*pdf_q1(n_pdfrw1))/
pdf(pdf_pdg_code1(n_pdfrw1),pdf_x1(n_pdfrw1),pdf_q1(n_pdfrw1-1))*
# Note the central scale reweighting by scalefact above
# (if n_pdfrw1 = 1, need to reweight the initial state pdf scale)
# Also note that no scale should be larger than the last one (including
# scalefact) for that beam.
#
# Now the same thing for all pdfs in beam 2 (n_pdfrw2)
pdf(pdf_pdg_code2(1),pdf_x2(1),pdf_q2(1))*       # initial state pdf
 ...
pdf(pdf_pdg_code2(n_pdfrw2),pdf_x2(n_pdfrw2),scalefact*pdf_q2(n_pdfrw2))/
pdf(pdf_pdg_code2(n_pdfrw2),pdf_x2(n_pdfrw2),pdf_q2(n_pdfrw2-1))*
/ total_reweight_factor                          # corr. factor from MG run
* 0 if (QCUT > SMIN or SMAX > max(QCUT, SCOMP)), otherwise 1}}}