Generating a Stop Simplified Model in Madgraph 5

The start-up procedure is like for Maxim's madgraph tutorial:

% cd <DIR>/MadGraph5_v1_4_8_2
% ./bin/mg5

Now, let's implement the simplified model within MadGraph 5's "mssm" model:

mg5>import model mssm
INFO: load particles 
INFO: load vertices 
INFO: Restrict model mssm with file models/mssm/restrict_default.dat . 
INFO: Run "set stdout_level DEBUG" before import for more information. 
INFO: Change particles name to pass to MG5 convention 
Kept definitions of multiparticles l- / j / vl / l+ / p / vl~ unchanged

Note that now, MG5 knows about many more particles than before!

mg5>display particles
Current model contains 48 particles:
w+/w- x1+/x1- x2+/x2- h+/h- ve/ve~ vm/vm~ vt/vt~ e-/e+ mu-/mu+ ta-/ta+ u/u~ c/c~ t/t~ d/d~ s/s~ b/b~ sve/sve~ svm/svm~ svt/svt~ el-/el+ mul-/mul+ ta1-/ta1+ er-/er+ mur-/mur+ ta2-/ta2+ ul/ul~ cl/cl~ t1/t1~ ur/ur~ cr/cr~ t2/t2~ dl/dl~ sl/sl~ b1/b1~ dr/dr~ sr/sr~ b2/b2~ 
a z g n1 n2 n3 n4 go h1 h2 h3 

Let's try to generate the process in which the lightest stop is pair-produced:

mg5>generate p p > t1 t1~ 
mg5>output madevent stopA
mg5>open index.html
mg5>launch
[to save time, modify run_card.dat so that it only generates 2000 events]

Now let's check a few things:

	•	First, look at the diagrams from index.html … note that QED diagrams weren't included (which is what we want), because of the settings of MadGraph.  If we *did* want to include QED processes and QCD, we could have run 

		mg5>generate p p > t1 t1~ QED=2 QCD=2
		mg5>output madevent stopEWK
		mg5>open index.html

	•	Now, what about the stop decays?  Here is a sample from the LHE file:
		
<event>
 4   0  0.6559200E-03  0.4496188E+03  0.7957747E-01  0.1037442E+00
       21   -1    0    0  501  502  0.00000000000E+00  0.00000000000E+00  0.57410529230E+03  0.57410529230E+03  0.00000000000E+00 0.  1.
       21   -1    0    0  502  503  0.00000000000E+00  0.00000000000E+00 -0.35980485976E+03  0.35980485976E+03  0.00000000000E+00 0.  1.
  1000006    1    1    2  501    0  0.19917133879E+03 -0.52468272527E+02  0.38930938239E+02  0.45130108538E+03  0.39966850000E+03 0.  0.
 -1000006    1    1    2    0  503 -0.19917133879E+03  0.52468272527E+02  0.17536949430E+03  0.48260906668E+03  0.39966850000E+03 0.  0.
</event>

	 	The stop has not decayed at all!  If we run Pythia, it will decay the stops, according to the mixing angles, etc. specified in param_card.dat.  (If you haven't done so yet, you can run:
		mg5> install pythia-pgs 
		to set up pythia)

So, how do we generate only the decays we want?  One option is to specify the decay in MadGraph: 

mg5>generate p p > t1 t1~ QED=0, t1~ > t~ n1 , t1 > t n1
mg5>output madevent stopB
mg5>open index.html
mg5>launch

When you run MadEvent, make sure you edit param_card.dat so that there's kinematic phase space for the decay!
e.g. in my default param_card.dat
      1000006 3.996685e+02 #   msu3
      1000024 1.816965e+02 #   mch1
=> it will run!

Note this will take a bit longer than the last process.  When it's done, compare the MadGraph reported cross-section in stopB to your stopA run -- it's a lot smaller!  The reason is that here MadGraph is reporting a cross-section times (branching ratio)^2 for the SUSY model whose parameters were entered in param_card.dat.

Two suggested exercises (on your own time)
	1.	In a Mathematica notebook, write an analysis script to "mock up" an ATLAS or CMS search for stops.  Choose a few mass points at which to compare your mock-up's acceptance to what the experiments report in their papers.  Note that to get even decent agreement, you would need:
	⁃	a good parametrization of efficiency for lepton identification and (especially) b-tagging
	⁃	rather than using MadGraph's parton-level events directly (treating each parton as a jet), you should run the Pythia parton shower and hadronization, using a jet algorithm similar to the experiments (usually anti-kT for BSM analyses, with cone size specified in the paper)
	⁃	you may also want to take into account the effect of jet energy resolution on jet energies and missing energy (this can be done crudely by using PGS or Delphes).  In my experience, the improvement from using PGS over "hadron-level MC truth" is usually not dramatic, but it's good to run both, as this gives you some idea of whether your mock-up's performance depends sensitively on your model of the detector.  If it DOES depend sensitively (i.e. your MC-truth efficiency and your PGS efficiency are very different), it's a sign that you should be very skeptical of your findings, because NEITHER PGS nor MC-truth is a faithful model of a real detector!)
	2.	Try generating some of the more complex simplified models we discussed, or some of your own.