MC propagation

Each simulation stage should only contain information from the previous one to avoid using MC information in a simulation stage where it will not be available later on with real data from the experiment. Nevertheless in many cases the MC information is important and a way to get access to the MC information from each simulation stage is needed.

The basis of the introduced MC propagation mechanism in FairRoot is that each stored data object contains the information which data was used to create it. For example each MC point knows from which MC track it was coming from or each reconstructed track knows which hit was used to create the track. Technically the information is stored in an array of FairLink objects. A FairLink is a pair of a branch ID, indicating from which data-branch the data is coming from, and an index, which indicates the position in the TCLonesArray of the data-branch. A FairLink uniquely identifies a data-object and allows the retrieving of the data-object from the root-file where it is stored in.

Following the FairLinks allows to get access from each simulation stage back to each previous stage to extract (not only) the MC information.
More information on how to use the FairLink data can be found here

How to implement MC propagation

IMPORTANT With SVN revision 9240 the implementation of MC propagation has changed! A FairLink consists now of a string and an integer. The string is the branch name, where the data is stored in. The integer is the position in the TClonesArray as usual. Internally the FairRootManager translates the branch name into a unique ID which is stored in the file.

In addition the setting of links now checks if the data the link is pointing to is transient. If this is the case the link information of the previous stage is saved instead of the original link. In this way it is ensured that the chain of link is not interrupted by intermediated data which is not stored.

  1. Derive your data class from FairMultiLinkedData*
    • If your data class already derives from FairHit or FairMCPoint this is not necessary
    • Be careful Multiple inheritance from different classes causes problems in ROOT so better avoid it.
  2. Use one of the methods SetLink, SetLinks, AddLink or AddLinks
    • Use the methods SetLink/SetLinks, AddLink/AddLinks to save the FairLink to the previous data-object
    • SetLink(FairLink) overrides an existing link array with one FairLink
    • SetLinks(FairMultiLinkedData) overrides an existing array with an array of FairLinks
    • AddLink(FairLink) / AddLinks(FairMultiLinkedData) adds the given data at the end of the existing FairLink array.
  3. Compile / Check/ Upload
    • Compile your code
    • Run a simulation
    • Check if your data class contains a fLinks branch with entries for fType and fIndex
    • Upload your code to SVN
    • Thats it. You have successfully implemented MC propagation into your code.

How to use MC Propagation


If you want to use MC Propagation for your data classes you first have to follow the description on Howto implement MC Propagation . After this you have to initialize your new data classes in the InitBranches method of the two tasks: PndMCMatchCreatorTask and PndMCMatchLoaderTask. This steps have to be done once for each new data class.

The base class which contains all the link informations is PndMCMatch. It is filled by the two Tasks: PndMCMatchCreatorTask and PndMCMatchLoaderTask

PndMCMatchCreatorTask
PndMCMatchCreatorTask runs through all initialized branches and extracts the stored FairLinks into the PndMCMatch class. After this the link information is stored in an individual branch so that it can be accessed without touching the original data. If this task is called in addition with other simulation tasks, it should be called as the last one, otherwise it will not contain the link information of the tasks coming after it. If you do not want to store the link data in the same root file as your data you should call PndMCMatchCreatorTask in an individual macro with a new output file name. You can work with this new file without loading the data files.

PndMCMatchLoaderTask
PndMCMatchLoaderTask reads in the links from a link root file or from a data file which contains a link branch. PndMCMatchLoaderTask is automatically called if you use the PndMCMatchCreatorTask.

After you have loaded the link information you can use either the GetMCInfo or the GetMCInfoSingle method of PndMCMatch to access the wanted link information:

GetMCInfo
PndMCResult GetMCInfo(TString start, TString stop) needs two parameters: the start branch from where you want to have your links and the stop branch where the links should point to. So if you want to know which MC track belongs to your reconstructed tracks start would be "PndTrack" and stop "MCTrack". The output is of type PndMCResult. Basically this is an array of FairLinks. The FairLinks are ordered in the same way as the start data. So if you want to know which MC tracks belong to your reconstructed track 2 you call GetEntry(2). Be careful: this can be more than one! So the output is a FairLink array which you can iterate by calling GetLink(index). These links come with a weight which tells you how much influence the link data had to your output. The weighting of links is quite complicated and will be covered in an own wiki page.

PndMCEntry GetMCInfoSingle(FairLink aLink, TString stop)has the same functionality as GetMCInfo. The only difference is that you start with an individual FairLink to a specific data object and you get only the FairLinks for this object. Let us take the same example like before: reconstructed track 2 -> MC tracks. Then your code would look like GetMCInfoSingle(FairLink(kTrack, 2), kMCTrack);.

Example


example.cxx

void PndMCTestMomentumCompare::Exec(Option_t* opt)
{
   PndMCResult myResult = fMCMatch->GetMCInfo("PndTrack", "MCTrack");
   std::cout << myResult;
   for (int i = 0; i < myResult.GetNEntries(); i++){
      PndMCEntry myLinks = myResult.GetEntry(i);
      PndTrack* myTrack = (PndTrack*)fTrack->At(i);
      std::cout << "TrackMatch for Track " << i << std::endl;
      std::cout << "P: " << myTrack->GetParamFirst().GetSDMomentum().Mag() << std::endl;
      std::cout << "Belongs to: " << std::endl;
      for (int j = 0; j < myLinks.GetNLinks(); j++){
         if (myLinks.GetLink(j).GetType() == kMCTrack){
            std::cout << "MCTrack " << myLinks.GetLink(j).GetIndex() << std::endl;
            PndMCTrack* myMCTrack = (PndMCTrack*)fMCTrack->At(myLinks.GetLink(j).GetIndex());
            std::cout << "P: " << myMCTrack->GetMomentum().Mag() << " PID: " << myMCTrack->GetPdgCode() << std::endl;
            std::cout << "--------------------------------" << std::endl;
         }
      }
      std::cout << std::endl;
   }
}