Step 6: Monte Carlo Tuning Sample! (Fortran Version)

Super-K Code Exercise: Making a Monte Carlo Tuning Sample

This exercise will allow you to employ skills learned from the previous tutorials, and in fact is a task similar to the kind of thing one often needs to do in reality. The idea is to create a sample of MC events that one can compare to data. You will take a file of stopping muon data events, fit (reconstruct) them, and create kinematics files using the fit results in order to simulate a matched MC sample.

Monte Carlo Tuning

Detector simulation code usually has a number of parameters that can be tweaked in order to change the nature of the simulation. Some examples are: sizes or densities of detector elements, optical properties of detector materials, properties of the electronics. Sometimes these parameters are quantities that can be known precisely, either from the way the detector is constructed or else independently measured: for instance, in Super-K the sizes of the PMTs are well known, as is the density of the water. Other parameters may not be so well known, either because it’s hard to measure them directly or they might change with time. For instance, transparency of the water in Super-K can vary with its cleanliness. For the case when the parameters aren’t well known, one can “tune” them: i.e. determine parameter values in order that simulated events match real data events as well as possible.

A typical way of tuning the MC is to obtain a data sample, and simulate a MC sample which in principle should match it. Then, one looks at various key distributions for data and Monte Carlo (the nature of which depends on what parameters you’re looking at) and tweaks the parameters until the distributions match are well matched.

A Stopping Muon Example

For this example, you will take a sample of selected stopping muons in Super-K (real data events). You will then fit them using a simple fitter called muboy. The fitter outputs a vertex and direction for each muon. You will use this information to create a kinematics file that corresponds to the stopping muon sample, i.e. each real data event will have a corresponding event in the kinematics file for simulation. Then, run skdetsim to simulate the events. Finally, you will create plots for your data and MC distributions for comparison. The point of the example is to create a sample that can be used for tuning; we’ll stop short of the actual tuning step (although some of you will be doing this later). I am not providing complete canned code for each of these steps: your challenge is to fill in the blanks to complete the task. Most of the examples we’ve done so far should give you guidance.

Files for this example

Copy the following files from

$TUTEXAMPLES/tutorials

into your working directory:

  • muboy_fit.F: contains sample source code for fitting muons.
  • The shell script for running the program, muboy_fit.sh.
  • $TUTEXAMPLES/patmue.run069486.002.zbs contains primarily stopping muon data events (superscan them to have a look). You need not copy this file into your working directory (it is a huge file). This is a data file from SK-IV.
  • sk4_odtune.card contains parameters for skdetsim appropriate to Super-K IV.

Basic Steps

  • First, examine muboy_fit.F. You will see it is similar to the previous examples, but has some extra variable declarations for muboy fit variables. For each event, it calls the muboy fitting routine, 
                   call muboy_zbs(nevsk,muboy_class,muboy_entry,
     +              muboy_dir,muboy_tracklength,
     +              muboy_goodness,muboy_numtracks,
     +              muboy_otherentry,muboy_numleft)

    All but the first of these arguments are the fitting results returned to you. Of these, variables you will need are:

    • muboy_class: this is the type of event. 1 means a through-going muon, and 2 means a stopper. The events in this example’s file will be overwhelmingly class 2.
    • muboy_entry: this is the entry point of the muon on the ID.
    • muboy_dir: the fit direction vector of the muon.
    • muboy_tracklength: the tracklength of the muon in the detector.

     

  • Compile and run the fitting program on the stopping muon file, and check the output written to the screen.
  • Now, modify muboy_fit.F to open a kinematics file and write an entry for each data event, where the vertex and direction of the event correspond to the vertex and direction resulting from the muboy fit. As a first pass, you can use the muboy entry point (in the ID) as the vertex. For the energy, take the pathlength provided by muboy and multiply by the energy loss of a muon in water, which is approximately 2 MeV/cm. You will probably need to consult a Fortran manual or tutorial to figure out how to write output in the proper format (alternatively, you can write to a simple unformatted text file, and use Perl or another language to write the properly formatted output, if you are more comfortable with that).
  • Once you have your kinematics file, simulate the events with skdetsim. You can check with superscan that the data and MC events correspond to each other.  Note that you will need to modify the maximum number of events simulated (VECT-NEVT) in the card file, the default is 25, you will likely want more.
  • Finally, write a program (I am not providing one– use the examples so far as a guide) to read the data file, read the MC file, and output the number of hits in the ID, nqisk, for comparison between data and MC.
  • Make a plot of the distribution of number of ID hits for data, with a plot of MC superimposed, using ROOT.
  • Note that there are parameters in the card file and in muboy_fit.F for the number of events to loop over.  If your data and MC don’t match up, check that you have the same number of events in data and MC.

A Refinement

Since the entry point given by muboy is on the ID cylinder, if you make the vertex the same as the muboy entry point, the muon will not be properly tracked in the OD (check this yourself with superscan). To simulate events with appropriate OD light, you will need to extrapolate the muboy track back to where it intersects the OD outer wall, and make the vertex the entry point in the OD.

Try this next, and make a plot of the distribution of OD hits for data and MC. Here are some tools that will help:

  • The routine crstnksk (standard SK software) finds the intersection of a line and a cylinder. You can find this in $SKOFL_ROOT/src/sklib/crstnksk.F: the documentation at the top of this source file tells you how to use it.
  • The include file geotnk.h has parameters which give the size of the SK tank ID and OD. You can find this in $SKOFL_ROOT/inc/geotnk.h.