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 /*
  17 $Log$
  18 Revision 1.4  2002/10/22 15:02:15  alibrary
  19 Introducing Riostream.h
  20
  21 Revision 1.3  2002/10/14 14:57:33  hristov
  22 Merging the VirtualMC branch to the main development branch (HEAD)
  23
  24 Revision 1.2.2.1  2002/10/14 09:45:57  hristov
  25 Updating VirtualMC to v3-09-02
  26
  27 Revision 1.2  2002/10/01 08:47:04  jchudoba
  28 Change loop order to run faster.
  29
  30 Revision 1.1  2002/09/17 08:37:12  jchudoba
  31 Classes to create and store tracks maps - correcpondence between track label and entry number in the TreeH
  32
  33 */
  34
  35 ////////////////////////////////////////////////////////////////////////
  36 //
  37 // AliTrackMapper.cxx
  38 // description:
  39 //        create an AliTrackMap - a relation between track label and
  40 //        it's index in TreeH. Check all branches with hits.
  41 //        Output is in the root file.
  42 //
  43 ////////////////////////////////////////////////////////////////////////
  44
  45 #include <Riostream.h>
  46
  47 #include "TTree.h"
  48 #include "TROOT.h"
  49 #include "TFile.h"
  50 #include "TBenchmark.h"
  51 #include "TStopwatch.h"
  52
  53 #include "AliDetector.h"
  54 #include "AliTrackMapper.h"
  55 #include "AliTrackMap.h"
  56 #include "AliRun.h"
  57 #include "AliHit.h"
  58
  59 ClassImp(AliTrackMapper)
  60
  61 //_______________________________________________________________________
  62 AliTrackMapper::AliTrackMapper():
  63   fDEBUG(0)
  64 {
  65   //
  66   // Default ctor
  67   //
  68 }
  69
  70
  71 //_______________________________________________________________________
  72 void AliTrackMapper::CreateMap(Int_t nEvents, Int_t firstEventNr,
  73                                const char* fnMap, const char* fnHits)
  74 {
  75   //
  76   // method to create a track map for a given number of events starting
  77   // from event firstEventNr. This method just opens and closes files and
  78   // loops over events, the creation of map is delegated to
  79   // CreateMap(Int_t eventNr, TFile* fileMap) method
  80   //
  81   TStopwatch timer;
  82   timer.Start();
  83
  84   TFile *fileMap=TFile::Open(fnMap,"new");
  85   if (!fileMap->IsOpen()) {cerr<<"Can't open output file "<<fnMap<<"!\n"; return;}
  86
  87   TFile *fileHits=TFile::Open(fnHits);
  88   if (!fileHits->IsOpen()) {cerr<<"Can't open input file "<<fnHits<<"!\n"; return;}
  89   if (!(gAlice=dynamic_cast<AliRun*>(fileHits->Get("gAlice")))) {
  90     cerr<<"gAlice have not been found on galice.root !\n";
  91     return;
  92   }
  93
  94   for (Int_t eventNr = firstEventNr; eventNr < firstEventNr+nEvents;
  95        eventNr++) {
  96     CreateMap(eventNr,fileMap);
  97   } // end loop over events
  98
  99   delete gAlice;
 100   gAlice = 0;
 101   fileHits->Close();
 102   delete fileHits;
 103   fileMap->Close();
 104   delete fileMap;
 105   timer.Stop();
 106   if (fDEBUG > 0) timer.Print();
 107 }
 108
 109 //_______________________________________________________________________
 110 Int_t  AliTrackMapper::CreateMap(Int_t eventNr, TFile* fileMap)
 111 {
 112   //
 113   // create an AliTrackMap for a given event
 114   // correct gAlice must be already present in memory
 115   //
 116   Int_t nGenPrimPlusSecParticles = gAlice->GetEvent(eventNr);
 117   if (fDEBUG > 1) cout<<"nGenPrimPlusSecParticles = "<<nGenPrimPlusSecParticles<<endl;
 118   if (nGenPrimPlusSecParticles < 1) {
 119     cerr<<"No primary particles found in event "<<eventNr<<endl;
 120     return -1;
 121   }
 122
 123   TTree *treeK = gAlice->TreeK();
 124   if (!treeK) {
 125     cerr<<"Error: Event "<<eventNr<<", treeK not found."<<endl;
 126     return -1;
 127   }
 128   Int_t nAllParticles = static_cast<Int_t>(treeK->GetEntries());
 129   Int_t *trackMap = new Int_t[nAllParticles];
 130   for (Int_t i = 0; i<nAllParticles; i++) {trackMap[i] = kNoEntry;}
 131
 132   TTree *treeH = gAlice->TreeH();
 133   if (!treeH) {
 134     cerr<<"Error: Event "<<eventNr<<", treeH not found."<<endl;
 135     return -1;
 136   }
 137   Int_t treeHEntries = static_cast<Int_t>(treeH->GetEntries());
 138   if (fDEBUG > 1) cout<<"treeHEntries "<<treeHEntries<<endl;
 139
 140   TObjArray *modules = gAlice->Detectors();
 141   if (!modules) {
 142     cerr<<"TObjArray with modules not found."<<endl;
 143     return -1;
 144   }
 145   Int_t nModules = static_cast<Int_t>(modules->GetEntries());
 146   AliHit* hit;
 147   for (Int_t treeHIndex = 0; treeHIndex < treeHEntries; treeHIndex++) {
 148     gAlice->ResetHits();
 149     treeH->GetEvent(treeHIndex);
 150     for (Int_t iModule = 0; iModule < nModules; iModule++) {
 151       AliDetector * detector = dynamic_cast<AliDetector*>(modules->At(iModule));
 152       if (!detector) continue;
 153       // process only detectors with shunt = 0
 154       if (detector->GetIshunt()) continue;
 155
 156       hit=dynamic_cast<AliHit*>(detector->FirstHit(-1));
 157       Int_t lastLabel=-1, label;
 158       for( ; hit; hit=dynamic_cast<AliHit*>(detector->NextHit()) ) {
 159         label=hit->Track();
 160         if (lastLabel != label) {
 161           if (label < 0 || label >= nAllParticles) {
 162             cerr<<"Error: label out of range. ";
 163             cerr<<"Event "<<eventNr<<" treeHIndex "<<treeHIndex<<" label = "<<label<<endl;
 164             return -2;
 165           }
 166           if (trackMap[label] >=0 && trackMap[label] != treeHIndex) {
 167             cerr<<"Error: different treeHIndex for label "<<label
 168                 <<" indeces: "<<trackMap[label]<<" != "<<treeHIndex;
 169             cerr<<" event "<<eventNr<<" detector "<<detector->GetName()<<endl;
 170             return -3;
 171           }
 172           trackMap[label] = treeHIndex;
 173           if (fDEBUG > 2) cout<<"treeHIndex, label = "<<treeHIndex<<" "<<label<<endl;
 174           lastLabel = label;
 175         }
 176       }
 177     }
 178   }
 179
 180   if (fDEBUG > 2) {
 181     for (Int_t i = 0; i < nAllParticles; i++) {
 182       cout<<eventNr<<"\t"<<i<<"\t"<<trackMap[i]<<endl;
 183     }
 184   }
 185   fileMap->cd();
 186   AliTrackMap* trackMapObject = new AliTrackMap(nAllParticles, trackMap);
 187   trackMapObject->SetEventNr(eventNr);
 188   trackMapObject->Write();
 189   delete trackMapObject;
 190
 191   delete [] trackMap;
 192   return 0;
 193 }
 194
 195 //_______________________________________________________________________
 196 AliTrackMap* AliTrackMapper::LoadTrackMap(Int_t eventNr, const char* fnMap, TFile* &fileMap) {
 197   //
 198   // read an AliTrackMap object for the given event eventNr from
 199   // the file fileMap
 200   //
 201   fileMap=TFile::Open(fnMap);
 202   if (!fileMap->IsOpen()) {cerr<<"Can't open file "<<fnMap<<" with map!\n"; return 0;}
 203   char mapName[20];
 204   sprintf(mapName,"AliTrackMap_%5.5d",eventNr);
 205   AliTrackMap* trackMapObject = dynamic_cast<AliTrackMap*>(fileMap->Get(mapName));
 206   if (!trackMapObject) {
 207     cerr<<"Error: map named "<<mapName<<" not found."<<endl;
 208     return 0;
 209   }
 210   return trackMapObject;
 211 }
 212
 213 //_______________________________________________________________________
 214 void AliTrackMapper::CheckTrackMap(Int_t eventNr, const char* fnMap) {
 215   //
 216   //
 217   //
 218   TFile *fileMap;
 219   AliTrackMap* trackMapObject = LoadTrackMap(eventNr, fnMap, fileMap);
 220   if (!trackMapObject) return;
 221
 222   trackMapObject->PrintValues();
 223
 224   delete trackMapObject;
 225   fileMap->Close();
 226   delete fileMap;
 227 }
 228