TPC/AliTPCTracking.C

   1 ////////////////////////////////////////////////////////////////////////
   2 //
   3 // AliTPCTracking.C
   4 //
   5 // date: 22.08.2002
   6 // author: Jiri Chudoba based on code of Jourij Belikov
   7 // version: 1.0
   8 // description:
   9 //      reconstructs of tracks in TPC inthe following steps:
  10 //         TPC cluster finding
  11 //         TPC track finding
  12 // input parameters:
  13 //        Int_t nEvents      ... nr of events to process
  14 //        Int_t firstEventNr ... first event number (starts from 0)
  15 //        const char* fileNameHits ... name of file with hits
  16 //        const char* fileNameDigits .. name of file with TPC digits
  17 //        const char* fileNameClusters .. name of file with TPC clusters (output)
  18 //        const char* fileNameTracks .. name of file with TPC tracks (output)
  19 //
  20 //        default file names correspond to pp production (2002-04)
  21 //
  22 // History:
  23 //
  24 //     18.03.2003 ... Char_t* replaced by const char*
  25 //
  26 //     03.03.2003 ... SetFieldFactor moved to AliTracker class and
  27 //                    LoadTPCParam moved to AliTPC class
  28 //                    TString replaced by Char_t*
  29 //
  30 //     20.11.2002 ... Changes due to a changed interface of AliTPCtracker.
  31 //                    Use Riostream.h instead of iostream.h
  32 //
  33 //     22.08.2002 ... first version
  34 //
  35 ////////////////////////////////////////////////////////////////////////
  36
  37 #if !defined(__CINT__) || defined(__MAKECINT__)
  38 #include "Riostream.h"
  39 #include "TFile.h"
  40 #include "TTree.h"
  41 #include "TBenchmark.h"
  42 #include "AliTPCtracker.h"
  43 #include "AliTPCtrackerMI.h"
  44 #include "AliTPCclusterer.h"
  45 #include "AliTPCclustererMI.h"
  46 #include "AliTPC.h"
  47 #include "AliRun.h"
  48 #include "AliHeader.h"
  49 #include "AliGenEventHeader.h"
  50 #include "AliTracker.h"
  51
  52 #endif
  53
  54 Int_t gDEBUG = 2;
  55
  56 Int_t AliTPCTracking(Int_t nEvents=1, Int_t firstEvent=0,
  57                      const char* fileNameHits="galice.root",
  58                      const char* fileNameDigits="tpc.digits.root",
  59                      const char* fileNameClusters="tpc.clusters.root",
  60                      const char* fileNameTracks="tpc.tracks.root",
  61                      Bool_t versionMI = kFALSE);
  62
  63 Int_t TPCFindClusters(Int_t nEvents=1, Int_t firstEvent=0,
  64                       const char* fileNameDigits="tpc.digits.root",
  65                       const char* fileNameClusters="tpc.clusters.root",
  66                       Bool_t versionMI = kFALSE);
  67 Int_t TPCFindClusters(Int_t nEvents, Int_t firstEvent,
  68                       TFile* fileDigits, TFile* fileClusters,
  69                       AliTPCParam* paramTPC=0);
  70 Int_t TPCFindClustersMI(Int_t nEvents, Int_t firstEvent,
  71                         TFile* fileDigits, TFile* fileClusters,
  72                         AliTPCParam* paramTPC=0);
  73 Int_t TPCFindTracks(Int_t nEvents=1, Int_t firstEvent=0,
  74                     const char* fileNameClusters="tpc.clusters.root",
  75                     const char* fileNameTracks="tpc.tracks.root",
  76                     Bool_t versionMI = kFALSE);
  77 Int_t TPCFindTracks(Int_t nEvents, Int_t firstEvent,
  78                     TFile* fileClusters, TFile* fileTracks,
  79                     AliTPCParam* paramTPC=0);
  80 Int_t TPCFindTracksMI(Int_t nEvents, Int_t firstEvent,
  81                       TFile* fileClusters, TFile* fileTracks,
  82                       AliTPCParam* paramTPC=0);
  83
  84 void FindVertex(Int_t iEvent, Double_t *vertex);
  85 void PrintVertex(TArrayF &primaryVertex);
  86
  87 ////////////////////////////////////////////////////////////////////////
  88 Int_t AliTPCTracking( Int_t nEvents, Int_t firstEvent,
  89                       const char* fileNameHits,
  90                       const char* fileNameDigits,
  91                       const char* fileNameClusters,
  92                       const char* fileNameTracks,
  93                       Bool_t versionMI) {
  94
  95   AliTracker::SetFieldFactor(fileNameHits,kFALSE);
  96
  97 // ********** Find TPC clusters *********** //
  98   if (fileNameDigits && fileNameClusters) {
  99     if (TPCFindClusters(nEvents,firstEvent,fileNameDigits,fileNameClusters,versionMI)) {
 100       cerr<<"Failed to get TPC clusters: !\n";
 101       return 1;
 102     }
 103   }
 104
 105 // ********** Find TPC tracks *********** //
 106   if (fileNameClusters && fileNameTracks) {
 107     if (TPCFindTracks(nEvents,firstEvent,fileNameClusters,fileNameTracks,versionMI)) {
 108       cerr<<"Failed to get TPC tracks !\n";
 109       return 2;
 110     }
 111   }
 112
 113   return 0;
 114 }
 115
 116 ////////////////////////////////////////////////////////////////////////
 117 Int_t TPCFindClusters(Int_t nEvents, Int_t firstEvent,
 118                       const char* fileNameDigits, const char* fileNameClusters,
 119                       Bool_t versionMI) {
 120
 121   Int_t rc;
 122   const Char_t *name="TPCFindClusters";
 123   if (gDEBUG>1) cout<<name<<" starts...\n";
 124   if (gDEBUG>1) gBenchmark->Start(name);
 125   TFile *fileClusters = TFile::Open(fileNameClusters,"recreate");
 126   TFile *fileDigits = TFile::Open(fileNameDigits);
 127   if (!fileDigits->IsOpen()) {
 128     cerr<<"Cannnot open "<<fileNameDigits<<" !\n";
 129     return 1;
 130   }
 131   if (!fileClusters->IsOpen()) {
 132     cerr<<"Cannnot open "<<fileNameClusters<<" !\n";
 133     return 1;
 134   }
 135
 136   if (versionMI) {
 137     rc = TPCFindClustersMI(nEvents,firstEvent,fileDigits,fileClusters);
 138   } else {
 139     rc = TPCFindClusters(nEvents,firstEvent,fileDigits,fileClusters);
 140   }
 141
 142   fileDigits->Close();
 143   fileClusters->Close();
 144   delete fileDigits;
 145   delete fileClusters;
 146   if (gDEBUG>1) gBenchmark->Stop(name);
 147   if (gDEBUG>1) gBenchmark->Show(name);
 148
 149   return rc;
 150 }
 151 ////////////////////////////////////////////////////////////////////////
 152 Int_t TPCFindClusters(Int_t nEvents, Int_t firstEvent,
 153                       TFile* fileDigits, TFile* fileClusters,
 154                       AliTPCParam* paramTPC) {
 155
 156   fileDigits->cd();
 157   if (!paramTPC) paramTPC = AliTPC::LoadTPCParam(fileDigits);
 158   if (!paramTPC) return 1;
 159
 160   for (Int_t iEvent = firstEvent; iEvent < firstEvent+nEvents; iEvent++){
 161     if (gDEBUG > 2) cout<<"TPCFindClusters: event "<<iEvent<<endl;
 162     AliTPCclusterer::Digits2Clusters(paramTPC, fileClusters, iEvent);
 163   }
 164   return 0;
 165 }
 166 ////////////////////////////////////////////////////////////////////////
 167 Int_t TPCFindClustersMI(Int_t nEvents, Int_t firstEvent,
 168                         TFile* fileDigits, TFile* fileClusters,
 169                         AliTPCParam* paramTPC) {
 170
 171   fileDigits->cd();
 172   if (!paramTPC) paramTPC = AliTPC::LoadTPCParam(fileDigits);
 173   if (!paramTPC) return 1;
 174
 175   AliTPCclustererMI clusterer;
 176   for (Int_t iEvent = firstEvent; iEvent < firstEvent+nEvents; iEvent++){
 177     if (gDEBUG > 2) cout<<"TPCFindClustersMI: event "<<iEvent<<endl;
 178     char treeName[100];
 179     sprintf(treeName, "TreeD_75x40_100x60_150x60_%d", iEvent);
 180     TTree* input = (TTree*) fileDigits->Get(treeName);
 181     fileClusters->cd();
 182     sprintf(treeName, "TreeC_TPC_%d", iEvent);
 183     TTree* output = new TTree(treeName, treeName);
 184     clusterer.SetInput(input);
 185     clusterer.SetOutput(output);
 186     clusterer.Digits2Clusters(paramTPC, iEvent);
 187   }
 188
 189   return 0;
 190 }
 191 ////////////////////////////////////////////////////////////////////////
 192 Int_t TPCFindTracks(Int_t nEvents, Int_t firstEvent,
 193                     const char* fileNameClusters, const char* fileNameTracks,
 194                     Bool_t versionMI) {
 195
 196   Int_t rc = 0;
 197   const Char_t *name="TPCFindTracks";
 198   if (gDEBUG>1) cout<<name<<" starts"<<endl;
 199   if (gDEBUG>1) gBenchmark->Start(name);
 200   TFile *fileTracks = TFile::Open(fileNameTracks,"recreate");
 201   TFile *fileClusters =TFile::Open(fileNameClusters);
 202
 203   if (versionMI) {
 204     rc = TPCFindTracksMI(nEvents, firstEvent, fileClusters, fileTracks);
 205   } else {
 206     rc = TPCFindTracks(nEvents, firstEvent, fileClusters, fileTracks);
 207   }
 208
 209   fileClusters->Close();
 210   fileTracks->Close();
 211   delete fileClusters;
 212   delete fileTracks;
 213   if (gDEBUG>1) gBenchmark->Stop(name);
 214   if (gDEBUG>1) gBenchmark->Show(name);
 215   return rc;
 216
 217 }
 218 ////////////////////////////////////////////////////////////////////////
 219 Int_t TPCFindTracks(Int_t nEvents, Int_t firstEvent,
 220                     TFile *fileClusters, TFile * fileTracks,
 221                     AliTPCParam* paramTPC) {
 222
 223   Int_t rc = 0;
 224   if (!paramTPC) paramTPC = AliTPC::LoadTPCParam(fileClusters);
 225   if (!paramTPC) return 1;
 226
 227   for (Int_t iEvent = firstEvent; iEvent < firstEvent+nEvents; iEvent++){
 228     if (gDEBUG > 2) cout<<"TPCFindTracks: event "<<iEvent<<endl;
 229     AliTPCtracker *tracker = new AliTPCtracker(paramTPC);
 230     tracker->SetEventNumber(iEvent);
 231     Double_t vertex[3];
 232     FindVertex(iEvent,vertex);
 233     tracker->SetVertex(vertex);
 234     fileClusters->cd();
 235     rc = tracker->Clusters2Tracks(0,fileTracks);
 236     delete tracker;
 237   }
 238   return rc;
 239 }
 240 ////////////////////////////////////////////////////////////////////////
 241 Int_t TPCFindTracksMI(Int_t nEvents, Int_t firstEvent,
 242                       TFile *fileClusters, TFile * fileTracks,
 243                       AliTPCParam* paramTPC) {
 244
 245   Int_t rc = 0;
 246   if (!paramTPC) paramTPC = AliTPC::LoadTPCParam(fileClusters);
 247   if (!paramTPC) return 1;
 248
 249   for (Int_t iEvent = firstEvent; iEvent < firstEvent+nEvents; iEvent++){
 250     if (gDEBUG > 2) cout<<"TPCFindTracksMI: event "<<iEvent<<endl;
 251     AliTPCtrackerMI tracker(paramTPC, iEvent);
 252     rc = tracker.Clusters2Tracks(0, fileTracks);
 253   }
 254   return rc;
 255 }
 256
 257 ////////////////////////////////////////////////////////////////////////
 258 void FindVertex(Int_t eventNr, Double_t *vertex) {
 259
 260   vertex[0] = vertex[1] = vertex[2] = 0.;
 261   if (!gAlice) {
 262     cerr<<"gAlice was not found! Using vertex position (0,0,0).\n";
 263     return;
 264   }
 265
 266   gAlice->GetEvent(eventNr);
 267   AliHeader *header = gAlice->GetHeader();
 268   if (!header) {
 269     cerr<<"header was not found!\n";
 270     return;
 271   }
 272   AliGenEventHeader* genEventHeader = header->GenEventHeader();
 273   if (!genEventHeader) {
 274     cerr<<"AliGenEventHeader was not found!\n";
 275     return;
 276   }
 277
 278   TArrayF primaryVertex(3);
 279   genEventHeader->PrimaryVertex(primaryVertex);
 280   PrintVertex(primaryVertex);
 281   vertex[0] = static_cast<Double_t>(primaryVertex[0]);
 282   vertex[1] = static_cast<Double_t>(primaryVertex[1]);
 283   vertex[2] = static_cast<Double_t>(primaryVertex[2]);
 284 //  delete header;
 285   delete genEventHeader;
 286   return;
 287
 288 }
 289 ////////////////////////////////////////////////////////////////////////
 290 void PrintVertex(TArrayF &primaryVertex)
 291 {
 292   cout <<"Vertex: "
 293        <<primaryVertex[0]<<" "
 294        <<primaryVertex[1]<<" "
 295        <<primaryVertex[2]<<" "<<endl;
 296 }
 297 ////////////////////////////////////////////////////////////////////////