STEER/AliTrackMapper.cxx

   1 /**************************************************************************
   2  * Copyright(c) 1998-1999, ALICE Experiment at CERN, All rights reserved. *
   3  *                                                                        *
   4  * Author: The ALICE Off-line Project.                                    *
   5  * Contributors are mentioned in the code where appropriate.              *
   6  *                                                                        *
   7  * Permission to use, copy, modify and distribute this software and its   *
   8  * documentation strictly for non-commercial purposes is hereby granted   *
   9  * without fee, provided that the above copyright notice appears in all   *
  10  * copies and that both the copyright notice and this permission notice   *
  11  * appear in the supporting documentation. The authors make no claims     *
  12  * about the suitability of this software for any purpose. It is          *
  13  * provided "as is" without express or implied warranty.                  *
  14  **************************************************************************/
  15
  16 /* $Header$ */
  17
  18 ////////////////////////////////////////////////////////////////////////
  19 //
  20 // AliTrackMapper.cxx
  21 // description:
  22 //        create an AliTrackMap - a relation between track label and
  23 //        it's index in TreeH. Check all branches with hits.
  24 //        Output is in the root file.
  25 //
  26 ////////////////////////////////////////////////////////////////////////
  27
  28 #include <Riostream.h>
  29
  30 #include "TTree.h"
  31 #include "TROOT.h"
  32 #include "TFile.h"
  33 #include "TBenchmark.h"
  34 #include "TStopwatch.h"
  35
  36 #include "AliDetector.h"
  37 #include "AliTrackMapper.h"
  38 #include "AliTrackMap.h"
  39 #include "AliRun.h"
  40 #include "AliHit.h"
  41
  42 ClassImp(AliTrackMapper)
  43
  44 //_______________________________________________________________________
  45 AliTrackMapper::AliTrackMapper():
  46   fDEBUG(0)
  47 {
  48   //
  49   // Default ctor
  50   //
  51 }
  52
  53
  54 //_______________________________________________________________________
  55 void AliTrackMapper::CreateMap(Int_t nEvents, Int_t firstEventNr,
  56                                const char* fnMap, const char* fnHits)
  57 {
  58   //
  59   // method to create a track map for a given number of events starting
  60   // from event firstEventNr. This method just opens and closes files and
  61   // loops over events, the creation of map is delegated to
  62   // CreateMap(Int_t eventNr, TFile* fileMap) method
  63   //
  64   TStopwatch timer;
  65   timer.Start();
  66
  67   TFile *fileMap=TFile::Open(fnMap,"new");
  68   if (!fileMap->IsOpen()) {cerr<<"Can't open output file "<<fnMap<<"!\n"; return;}
  69
  70   TFile *fileHits=TFile::Open(fnHits);
  71   if (!fileHits->IsOpen()) {cerr<<"Can't open input file "<<fnHits<<"!\n"; return;}
  72   if (!(gAlice=dynamic_cast<AliRun*>(fileHits->Get("gAlice")))) {
  73     cerr<<"gAlice have not been found on galice.root !\n";
  74     return;
  75   }
  76
  77   for (Int_t eventNr = firstEventNr; eventNr < firstEventNr+nEvents;
  78        eventNr++) {
  79     CreateMap(eventNr,fileMap);
  80   } // end loop over events
  81
  82   delete gAlice;
  83   gAlice = 0;
  84   fileHits->Close();
  85   delete fileHits;
  86   fileMap->Close();
  87   delete fileMap;
  88   timer.Stop();
  89   if (fDEBUG > 0) timer.Print();
  90 }
  91
  92 //_______________________________________________________________________
  93 Int_t  AliTrackMapper::CreateMap(Int_t eventNr, TFile* fileMap)
  94 {
  95   //
  96   // create an AliTrackMap for a given event
  97   // correct gAlice must be already present in memory
  98   //
  99   Int_t nGenPrimPlusSecParticles = gAlice->GetEvent(eventNr);
 100   if (fDEBUG > 1) cout<<"nGenPrimPlusSecParticles = "<<nGenPrimPlusSecParticles<<endl;
 101   if (nGenPrimPlusSecParticles < 1) {
 102     cerr<<"No primary particles found in event "<<eventNr<<endl;
 103     return -1;
 104   }
 105
 106   TTree *treeK = gAlice->TreeK();
 107   if (!treeK) {
 108     cerr<<"Error: Event "<<eventNr<<", treeK not found."<<endl;
 109     return -1;
 110   }
 111   Int_t nAllParticles = static_cast<Int_t>(treeK->GetEntries());
 112   Int_t *trackMap = new Int_t[nAllParticles];
 113   for (Int_t i = 0; i<nAllParticles; i++) {trackMap[i] = kNoEntry;}
 114
 115   TTree *treeH = gAlice->TreeH();
 116   if (!treeH) {
 117     cerr<<"Error: Event "<<eventNr<<", treeH not found."<<endl;
 118     return -1;
 119   }
 120   Int_t treeHEntries = static_cast<Int_t>(treeH->GetEntries());
 121   if (fDEBUG > 1) cout<<"treeHEntries "<<treeHEntries<<endl;
 122
 123   TObjArray *modules = gAlice->Detectors();
 124   if (!modules) {
 125     cerr<<"TObjArray with modules not found."<<endl;
 126     return -1;
 127   }
 128   Int_t nModules = static_cast<Int_t>(modules->GetEntries());
 129   AliHit* hit;
 130   for (Int_t treeHIndex = 0; treeHIndex < treeHEntries; treeHIndex++) {
 131     gAlice->ResetHits();
 132     treeH->GetEvent(treeHIndex);
 133     for (Int_t iModule = 0; iModule < nModules; iModule++) {
 134       AliDetector * detector = dynamic_cast<AliDetector*>(modules->At(iModule));
 135       if (!detector) continue;
 136       // process only detectors with shunt = 0
 137       if (detector->GetIshunt()) continue;
 138
 139       hit=dynamic_cast<AliHit*>(detector->FirstHit(-1));
 140       Int_t lastLabel=-1, label;
 141       for( ; hit; hit=dynamic_cast<AliHit*>(detector->NextHit()) ) {
 142         label=hit->Track();
 143         if (lastLabel != label) {
 144           if (label < 0 || label >= nAllParticles) {
 145             cerr<<"Error: label out of range. ";
 146             cerr<<"Event "<<eventNr<<" treeHIndex "<<treeHIndex<<" label = "<<label<<endl;
 147             return -2;
 148           }
 149           if (trackMap[label] >=0 && trackMap[label] != treeHIndex) {
 150             cerr<<"Error: different treeHIndex for label "<<label
 151                 <<" indeces: "<<trackMap[label]<<" != "<<treeHIndex;
 152             cerr<<" event "<<eventNr<<" detector "<<detector->GetName()<<endl;
 153             return -3;
 154           }
 155           trackMap[label] = treeHIndex;
 156           if (fDEBUG > 2) cout<<"treeHIndex, label = "<<treeHIndex<<" "<<label<<endl;
 157           lastLabel = label;
 158         }
 159       }
 160     }
 161   }
 162
 163   if (fDEBUG > 2) {
 164     for (Int_t i = 0; i < nAllParticles; i++) {
 165       cout<<eventNr<<"\t"<<i<<"\t"<<trackMap[i]<<endl;
 166     }
 167   }
 168   fileMap->cd();
 169   AliTrackMap* trackMapObject = new AliTrackMap(nAllParticles, trackMap);
 170   trackMapObject->SetEventNr(eventNr);
 171   trackMapObject->Write();
 172   delete trackMapObject;
 173
 174   delete [] trackMap;
 175   return 0;
 176 }
 177
 178 //_______________________________________________________________________
 179 AliTrackMap* AliTrackMapper::LoadTrackMap(Int_t eventNr, const char* fnMap, TFile* &fileMap) {
 180   //
 181   // read an AliTrackMap object for the given event eventNr from
 182   // the file fileMap
 183   //
 184   fileMap=TFile::Open(fnMap);
 185   if (!fileMap->IsOpen()) {cerr<<"Can't open file "<<fnMap<<" with map!\n"; return 0;}
 186   char mapName[20];
 187   sprintf(mapName,"AliTrackMap_%5.5d",eventNr);
 188   AliTrackMap* trackMapObject = dynamic_cast<AliTrackMap*>(fileMap->Get(mapName));
 189   if (!trackMapObject) {
 190     cerr<<"Error: map named "<<mapName<<" not found."<<endl;
 191     return 0;
 192   }
 193   return trackMapObject;
 194 }
 195
 196 //_______________________________________________________________________
 197 void AliTrackMapper::CheckTrackMap(Int_t eventNr, const char* fnMap) {
 198   //
 199   //
 200   //
 201   TFile *fileMap;
 202   AliTrackMap* trackMapObject = LoadTrackMap(eventNr, fnMap, fileMap);
 203   if (!trackMapObject) return;
 204
 205   trackMapObject->PrintValues();
 206
 207   delete trackMapObject;
 208   fileMap->Close();
 209   delete fileMap;
 210 }
 211