Begun: 2004/10/04, M. Kordosky kordosky@fnal.gov
Last significant change (new ion codes): 2007/07/16, R. Hatcher rhatcher@fnal.gov
Modified to add ParticleTransportSim section: 2007/04/22, S. Kasahara schubert@hep.umn.edu
This section documents some of the Monte Carlo truth information which may be accessed by offline jobs. The emphasis is on the StdHep block. The information was gathered from emails, GEANT documentation, private conversations and ``legacy'' web-pages (most of which were authored by Robert Hatcher).
It is possible to configure GMINOS to record a list of interesting processes in each event. The information eventually shows up as a TClonesArray of TParticle objects stored inside SimSnarlRecords. Inside the offline code (e.g., in a JobModule) users may retrieve the information via:
SimSnarlRecord* simrec = dynamic_cast<SimSnarlRecord*>
(mom->GetFragment("SimSnarlRecord"));
const TClonesArray* simstdheparray = dynamic_cast<const TClonesArray*>
(simrec -> FindComponent("TClonesArray",
"StdHep"));
// loop over array of TParticles
Int_t nstdhep = simstdheparray->GetLast()+1;
for ( Int_t istd = 0; istd < nstdhep; istd++ ) {
TParticle* simstdhep =
dynamic_cast<TParticle*>(simstdheparray->At(istd));
// Do whatever you want
// See NtpMCModule::FillNtpMCStdHep() for examples
}
In addition the truth information is recorded by the JobControl module
NtpMCModule as NtpMCStdHep objects which are stored
in NtpMCRecords. That's a mouthful! Exploring the
information in NtpMCRecords is pretty easy inside a loon session:
// open a file
loon [0] TFile f("f24110001_0000.sntp.R1.9.root");
// what's in the file?
loon [1] .ls
TFile** f24110001_0000.sntp.R1.9.root Persistency file
TFile* f24110001_0000.sntp.R1.9.root Persistency file
KEY: TTree NtpMC;1 Persistency stream tree
KEY: TTree NtpSR;1 Persistency stream tree
KEY: TTree NtpTH;1 Persistency stream tree
KEY: TNamed FileEnd;1 PerFile End Key
// what do NtpMC records look like?
loon [2] NtpMC->Show(0);
======> EVENT:0
NtpMCRecord = NULL
fUniqueID = 0
fBits = 50331648
fName =
fTitle =
fHeader.fUniqueID = 0
fHeader.fBits = 50331648
fHeader.fVldContext.fUniqueID = 0
fHeader.fVldContext.fBits = 50331648
fHeader.fVldContext.fDetector = 2
fHeader.fVldContext.fSimFlag = 4
fHeader.fVldContext.fTimeStamp.fSec = 1079981341
fHeader.fVldContext.fTimeStamp.fNanoSec = 51
fHeader.fRun = 24110001
fHeader.fSubRun = 1
fHeader.fRunType = 0
fHeader.fErrorCode = 0
fHeader.fSnarl = 0
fHeader.fTrigSrc = 0
fHeader.fTimeFrame = -1
fHeader.fEvent = -1
mc = 1
mc.fUniqueID = 0
mc.fBits = 50331648
mc.index = 0
mc.inu = 12
mc.inunoosc = 14
mc.itg = 2212
mc.iboson = -99999
mc.iresonance = 1003
mc.iaction = 0
mc.a = 56.000000
mc.z = 26.000000
mc.sigma = 216.959991
mc.x = 0.056173
mc.y = 0.337257
mc.q2 = -0.134168
mc.w2 = 3.070345
mc.emfrac = 0.000000
mc.vtxx = 169.134476
mc.vtxy = -101.748550
mc.vtxz = 959.161987
mc.p4neu[4] = 0.000005 , 0.220568 , 3.796557 , 3.802959
mc.p4neunoosc[4] = 0.000005 , 0.220566 , 3.796557 , 3.802959
mc.p4tgt[4] = 0.201376 , -0.069906 , 0.001095 , 0.928162
mc.p4shw[4] = 0.266551 , 0.203926 , 1.338785 , 1.321045
mc.p4mu1[4] = 0.000000 , 0.000000 , 0.000000 , 0.000000
mc.p4mu2[4] = 0.000000 , 0.000000 , 0.000000 , 0.000000
mc.p4el1[4] = 0.000000 , 0.000000 , 0.000000 , 0.000000
mc.p4el2[4] = 0.000000 , 0.000000 , 0.000000 , 0.000000
mc.p4tau[4] = 0.000000 , 0.000000 , 0.000000 , 0.000000
stdhep = 9
stdhep.fUniqueID = 0, 0, 0, 0, 0, 0, 0, 0, 0
stdhep.fBits = 50331648, 50331648, 50331648, 50331648,
50331648, 50331648, 50331648, 50331648,
50331648
stdhep.index = 0, 1, 2, 3, 4, 5, 6, 7, 8
stdhep.parent[2] = -1 , -1
, -1 , -1
, 1 , 1
, 1 , 1
, 0 , 0
, 2 , 2
, 5 , 5
, 5 , 5
, 5 , 5
stdhep.child[2] = 4 , 4
, 2 , 3
, 5 , 5
, -1 , -1
, -1 , -1
, 6 , 8
, -1 , -1
, -1 , -1
, -1 , -1
stdhep.IstHEP = 0, 0, 11, 1, 1, 3, 1, 1, 1
stdhep.IdHEP = 12, 1056026000, 2212, 1055025000, 12,
10000222, -211, 211, 2212
stdhep.mass = 0.000000, 52.103485, 0.903351, 51.174881,
0.000000, 1.752240, 0.139568, 0.139568,
0.938280
stdhep.p4[4] = 0.000005 , 0.220568 , 3.796557 , 3.802959
, 0.000000 , 0.000000 , 0.000000 , 52.103485
, 0.201376 , -0.069906 , 0.001095 , 0.928162
, -0.201376 , 0.069906 , -0.001095 , 51.175323
, -0.266546 , 0.016642 , 2.457773 , 2.472240
, 0.467927 , 0.134020 , 1.339880 , 2.258881
, 0.093360 , 0.191389 , 0.293592 , 0.388615
, 0.254709 , 0.230348 , 0.752490 , 0.838843
, 0.119858 , -0.287717 , 0.293798 , 1.031423
stdhep.vtx[4] = 1.691345 , -1.017485 , 9.591620 , 0.000000
, 1.691345 , -1.017485 , 9.591620 , 0.000000
, 1.691345 , -1.017485 , 9.591620 , 0.000000
, 1.691345 , -1.017485 , 9.591620 , 0.000000
, 1.691345 , -1.017485 , 9.591620 , 0.000000
, 1.691345 , -1.017485 , 9.591620 , 0.000000
, 1.691345 , -1.017485 , 9.591620 , 0.000000
, 1.691345 , -1.017485 , 9.591620 , 0.000000
, 1.691345 , -1.017485 , 9.591620 , 0.000000
The NtpMC tree contains NtpMCRecord objects. The objects have three main components:
Two codes are used to denote the sort of interaction that occurred between the neutrino and the target. These codes, known as ``Iaction'' and ``Iresonance'', are shown in Tables 21.1-21.2. The tables are valid for NEUGEN3 but perhaps not for earlier versions.
| ||||||||||||||
| ||||||||||||||||||||||||||
Table 21.3 lists those particles most commonly seen in our events. The numbering convention and the full table are documented in the PDG. In particular many resonances have been left out of the table. Anti-particles have a negative sign in front of the numeric code (e.g. a negative pion is -211). Geantinos are used by NEUGEN for bookkeeping.
Atomic nuclei or ions are marked with a 10 digit code, such as 1056026000. There are now two canonical encoding schemes. The older StdHep scheme uses 1AAAZZZJJJ where AAA is the atomic mass, ZZZ is the charge and JJJ is the (angular) spin. Currently the spin portion is unused. Therefore 1056026000 denotes Fe-56. The newer scheme was introduced in the 2006 PDG section on Monte Carlo Numbering (and hereforth referred to as PDG2006). This scheme uses 10LZZZAAAI where AAA and ZZZ are as above, L is the number of charmed Lambdas, and I distinguishes various excited states. Since MINOS is working with non-exotic matter both L and I will generally be zero. In the PDG2006 scheme Fe-56 would be 1000260560.
To facilitate working with these codes there are some useful functions defined in Util/UtilPDG.h; these are within the UtilPDG namespace (and thus need to be qualified).
The IstHEP status codes coming out of the NEUGEN3 generator deviate from the STDHEP conventions and are shown in Table 21.4. Particles of status code 1 are those that are fed into GEANT for propogation through the detector; these and their secondaries give rise to hits that record energy depositions in the active detector elements. The special status code 999 indicates a bogus entry used by GMINOS to transmit information about the flux; do not use these in any kinematics calculations.
In addition, if GMINOS deems the results of a process step to be ``interesting'' it will record the daughter products with a status code in the ranges of 201-230, 301-309, 1201-1230 where the code indicates the process that generated the secondaries as shown in Table 21.5. In general this GEANT3 tagging means a secondary exceeded a kinematic energy threshold set for that process (see Section 21.3). If the secondaries arise from a decay then the state of the decaying particle is also recorded. For decays (other than taus and heavy hadrons) the resulting daughters will be recorded as status code 205; in addition the parent that decayed will have an entry made with status code 1205 representing the state of the parent at the time of the decay (starting with daikon release).
| ||||||||||||||
| |||||||||||||||||||||||||||||||||||||||||||||||||||||||||||||||||||||||||||||||||||||||||||||||||||||||||||||||||||||||||||
Figures 21.1-21.4 exemplify various features that are recorded in the StdHep structure. These were generated using a JobPath that included RerootToTruth::Ana. The first column is the index into the array of TParticles, the second (in parentheses) is the status code. The ``parents'' indices delimit a range of indices that are the source of the particle. In the offline, using C++, indices start at 0 and thus -1 is used to indicate no linkage. The range of children (if any) are shown in square brackets. If one is looking for the the ``final state'' particles that come from the generator one should look for status code 1 and not follow daughter links. For each entry the right hand side of the first line shows the 4-momentum, and the second shows the position and time.
Figure 21.1 shows at 
CC event. The 
is a prompt
decay, so in MINOS there is no realistic possibility of seeing the tau
as an actual track. This is handled in GMINOS by allowing the tau
to propogate (starting position is entry 4, status code 2 in this dump),
recording the location of the decay (entry 6, status code 2) and then
putting the daughter particles in as status code 1 (entries 7-9). In this
event several particles underwent internuclear scattering (status code 14).
Figure 21.2 shows an ``opposite-sign di-muon'' event where the
charged current event contains a charmed particle that decays into a
muon (always of opposite charge as the primary). Like the CC
event the decay of interest is prompt leading to two entries for the
charmed particle (D0 meson in this case, entries 6 and 18, with status
code 2) and the decay products entered with status code 1.
Figure 21.3 shows the decay in flight of a single muon. The status code 205 indicates that this is em extra information beyond the primary event as desribed in Section 21.2; 205 in particular means a decay as documented in Table 21.5. Using these entries in kinematic calculations would be double or triple counting The first entry in each dump was the particle passed to GMINOS, while the second records the state of the muon at the point of the decay. For muons that decay in flight there is no threshold for the kinetic energy of the daughters; in the first case, the muon has essentially come to rest due to energy loss before the decay.
The final event, Figure 21.4, again demonstrates the GEANT
status codes. In this case one of the 
's decayed while being propogated
throught the detector with sufficient kinematic energy of a decay product
to warrent recording.
SAVT 1 7 28
would output the StdHep block for events one, seven and twenty-eight. One can save a range of events by using negative numbers:
SAVT -1 -28
That line would output the StdHep block for events one through twenty-eight. To write out every event you might do:
SAVT -1 -9999999
The C2ND card is used to set a momentum threshold (in GeV/c) for each GEANT process (listed in descending order from 1 to 39 in Table 21.5). Positive thresholds values will save all secondaries produced via the specified process provided that at least one secondary is above the threshold. Negative thresholds will save only those secondaries with a momentum larger than the (absolute value of) the threshold value. For example:
C2ND 12=0.01 C2ND 13=-0.05would save all secondaries from 'HADR' (code 12) processes in which one secondary had a momentum larger than 10MeV/c but for 'ECOH' (code 13) processes would only save secondaries with a momentum larger than 50MeV/c. A range of values can be given as:
C2ND 60*0.02 C2ND 12=-0.1which initially sets all sixty thresholds to 20MeV/c and then resets the 'HADR' threshold to 100MeV/c.
Currently the C2ND card only saves secondaries but in the future it may be helpful to save the primaries as well.
GEANT operates by propagating individual particles through the detector geometry - in other words it makes ``tracks''. Interactions may create new tracks. GEANT keeps track of the tracks (clever eh?) by assigning a unique number, called TrackId, to each track. The tracks do not have to appear in the StdHep output in order to be tagged with a TrackId. DigiScintHits record the ID of the track that created them and so propagate the information into the offline framework where it is used by, among other things, the Truthifier code. At one point in the past the author noticed negative TrackIds which made him fear that there was a bug. As it turns out, it's a feature, as Robert explains:
TrackId represents an index into the StdHep table - but that table would grow to outrageous size if every particle were recorded. So instead we use the convention that -TrackId means a hit from a particle that the true TrackId begat (at some level above). Thus deep in an EM shower initiated by a high energy electron (TrackId=5 for instance) some lowly electron is going to have TrackId=-5 but the state of all the e+/e-'s in the shower are not recorded.
More recently someone noticed that hits created about seven micro-seconds after the start of the event were being attributed, via TrackId, to a pion produced at the (neutrino) vertex. In light of Robert's comment above, it's likely that the hits were produced by the neutron progeny of the original pion.
Last significant change: 2007/04/22, S. Kasahara schubert@hep.umn.edu
ParticleTransportSim (PTSim) is the name for the package that handles the particle transport through the detector geometry. PTSim begins with the event initial state and results in the energy deposition hits and storage of selected secondaries generated during the simulation through the detector. It is the replacement for GMINOS particle transport beginning with the Eggplant Monte Carlo production.
The purpose of this section is to document conventions specific to PTSim.
The IstHEP status code assigned to secondaries by PTSim follows the convention of using UtilIstHEP::kProdMethodOffset (400) + the UtilIstHEP::EProdMethod code assigned at particle production. A listing of such process codes is shown in Table 21.5.1.
In addition, a special IstHEP status code of 1000 + UtilIstHEP::kProdMethodOffset + UtilIstHEP::kMDecay is reserved for the parent of a decay process, where the state of the parent has been recorded at the position of the decay. This is shown in Table 21.5.1.
The PTSimModule is configurable to allow the user to select the type of secondary to store to the output stdhep array. These configuration options, as well as the default behavior, are summarized in this section.
As with GMINOS, event interaction primaries, those with IstHEP status
code 
, are always stored to the output stdhep array.
The configuration parameter StdHepSave can be used to denote the snarl numbers for which selected secondary particles will be stored to the output StdHep array. The snarl numbers can be specified by range or as an itemized list. For example:
![\begin{figure}{\footnotesize\begin{verbatim}ihep stat type <parents > px py pz...
...1966
0 [ -1 -1] 1.56418 -5.88815 -14.022 0.0000e+00\end{verbatim}}
\end{figure}](img48.png)
![\begin{figure}{\footnotesize\begin{verbatim}ihep stat type <parents > px py pz...
...557
14 [ -1 -1] 1.00203 0.998222 1.08182 3.5489e-10\end{verbatim}}
\end{figure}](img50.png)
![\begin{figure}{\footnotesize\begin{verbatim}ihep stat type <parents > px py pz...
...3
111 [ -1 -1] -3.87914 -1.32384 21.6241 1.3886e-11\end{verbatim}}
\end{figure}](img51.png)
jc.Path("Reco").Mod("PTSimModule").Cmd("StdHepSave 1 4 10");
will store secondaries for snarls 1, 4, and 10.
Using a negative index will store the secondaries for a range of snarls. For example:
jc.Path("Reco").Mod("PTSimModule").Cmd("StdHepSave -1 -10");
stores secondaries for snarls 1 to 10, and:
jc.Path("Reco").Mod("PTSimModule").Cmd("StdHepSave 0 -5");
stores secondaries for snarls 0 to 5.
To store secondaries for every snarl, use:
jc.Path("Reco").Mod("PTSimModule").Cmd("StdHepSave 0 -9999999");
The default is to store secondaries for no snarls.
Note that if a secondary is selected for storage by one of the above means, the secondary's direct line of ancestors will be saved as well.
For example, to set both selection modes:
std::ostringstream oss;
Int_t selectmask = PTSim::kHit | PTSim::kMomentum;
oss << "StdHepSelectMask " << selectmask;
jc.Path("Reco").Mod("PTSimModule").Set(oss.str().c_str());
or to set selection by momentum threshold only one would specify:
std::ostringstream oss;
Int_t selectmask = PTSim::kMomentum;
oss << "StdHepSelectMask " << selectmask;
jc.Path("Reco").Mod("PTSimModule").Set(oss.str().c_str());
By default, StdHepSelectMask is set to PTSim::kMomentum.
The momentum threshold applied to secondaries can be adjusted using the StdHepThr parameter. For example, to specify a threshold of 0.20 GeV/c, use:
jc.Path("Reco").Mod("PTSimModule").Set("StdHepThr=0.20");
The default StdHepThr momentum threshold is 0.15 GeV/c.
Note that this threshold
is only used if the user has enabled selection mode PTSim::kMomentum
in the StdHepSelectMask.
The momentum threshold can also be adjusted by secondary ``type'' using the StdHepThrByType parameter. The ``type'' is the production mechanism and/or the PDG id of the secondary particle, where the production mechanism is one of the list of UtilIstHEP::EProdMethod enumerated mechanisms shown in Table 21.5.1.
The format for specifying a momentum threshold by type is:
jc.Path("Reco").Mod("PTSimModule")
.Cmd("StdHepThrByType <pthresh(GeV/c)> <mechanism type> <particle id>");
A momentum threshold specified by type will override the default momentum
threshold set by StdHepThr for the specified type.
For example, to specify momentum threshold 0.01 GeV/c for those particles produced by Bremsstrahlung:
std::ostringstream oss;
oss << "StdHepThrByType " << 0.01 << " " << UtilIstHEP::kMBrem;
jc.Path("Reco").Mod("PTSimModule").Cmd(oss.str().c_str());
As a second example, to specify momentum threshold 0.01 GeV/c for positrons produced by any mechanism, one would use:
std::ostringstream oss;
oss << "StdHepThrByType "<< 0.01 <<" "<< UtilIstHEP::kMNoProcess <<" "<< -11;
jc.Path("Reco").Mod("PTSimModule").Cmd(oss.str().c_str());
The third and final example specifies momentum threshold 0.30 for electrons produced by Delta Ray production:
std::ostringstream oss;
oss << "StdHepThrByType "<< 0.30 <<" "<< UtilIstHEP::kMDeltaRay <<" "<< 11;
jc.Path("Reco").Mod("PTSimModule").Cmd(oss.str().c_str());
By default, the StdHepThrByType momentum thresholds are set to zero for:
The parameter StdHepHitThr can be used to adjust the hit energy threshold for determining when a secondary resulting in an energy deposition hit is to be selected for storage. For example, to set the hit energy threshold to 0.002 GeV one would use:
jc.Path("Reco").Mod("PTSimModule").Set("StdHepHitThr=0.002");
The default hit threshold is 0.001 GeV.
Note that this threshold
is only used if the user has enabled selection mode PTSim::kHit
in the StdHepSelectMask.
The parameter StdHepSaveSibling can be used to specify if siblings of selected secondaries should also be selected for storage to the output StdHep array. Siblings are particles generated in the same production mechanism event, such as decay products from a single decay. For example, use:
jc.Path("Reco").Mod("PTSimModule").Set("StdHepSaveSibling=0");
to disable the storage of selected secondary sibilings. The default
value is 1 (enabled).