macros/CompareTime.C

   1 //
   2 // This macro draws time resolution for true hypothesis and fits it with gausian
   3 // on a given reference plane in a given momentum range.
   4 // The results is then printed to a PostScript file.
   5 // The filename is composed out of input parameters.
   6 //
   7 // Input parameters:
   8 // nRefPlane - indicate the reference plane.
   9 //   1: in TPC
  10 //   2: out TPC
  11 //   3: out TRD
  12 //
  13 // minP, maxP - total momentum cuts. momentum from gAlice are used for cuts
  14 //
  15 // The macro is using SORTED reference points in a file "trackRef.root".
  16 // Track references can be sorted using CopyReferences.C
  17 //
  18 // Before using check filenames and names of trees invlolved
  19 //
  20 // ITS tracks (nRefPlane == 1) are compared at inner TPC ref Plane
  21 // make shure they are properly propagated to this plane
  22 //
  23 // Sylwester Radomski, e-mail: S.Radomski@gsi.de
  24 // Feb 14, 2002
  25 //
  26
  27
  28 CompareTime (Int_t nRefPlane, Double_t minP, Double_t maxP)
  29 {
  30
  31   // Check input data Consistancy
  32   if (nRefPlane < 1 || nRefPlane > 3) {
  33     cout << "Wrong Reference Plane Id ( " << nRefPalne << ") [1-3] " << endl;
  34     return;
  35   }
  36
  37   if (maxP < minP) {
  38     cout << "MaxP lesser than MinP" << endl;
  39     return;
  40   }
  41
  42   // Style
  43
  44   gROOT->SetStyle("Plain");
  45   gStyle->SetOptStat(1110);
  46   gStyle->SetOptFit(11);
  47
  48   Double_t cc = 0.299792458;
  49   AliKalmanTrack::SetConvConst(100/cc/0.4);
  50
  51   const char* gFileName = "galice.root";
  52   const char* refFileName = "trackRef.root";
  53
  54   const char* trackFileName[] = {"AliTPCtracks.root", "AliTPCBackTracks.root", "AliTRDtracks.root"};
  55   const char *treeName[] = {"TreeT_ITSb_0", "seedsTPCtoTRD_0", "TRDb_0"};
  56   const char *refTreeName[] = {"TPC", "TPC", "TRD"};
  57
  58   // Histograms
  59
  60   const Int_t nTypes = 5;
  61   Int_t codes[] = {11, 13, 211, 321, 2212};
  62
  63   //Double_t cutoff[] = {20, 20, 20, 100, 100};
  64   Double_t cutoff[] = {30, 30, 30, 100, 200};
  65
  66   const char *names[nTypes] = {"Electrons", "Muons", "Pions", "Kaons", "Protons"};
  67
  68   TH1D *hDiffSame[nTypes];;
  69
  70   for(Int_t i=0; i<nTypes; i++)
  71     hDiffSame[i] = new TH1D(names[i], names[i], 50, -cutoff[i], cutoff[i]);
  72
  73   TH1D *hLength = new TH1D("length", "Length Resolution", 100, -0.5, 0.5);
  74
  75   // Particles
  76
  77   TFile *refFile = new TFile(gFileName, "READ");
  78   gAlice = (AliRun*)refFile->Get("gAlice");
  79   gAlice->GetEvent(0);
  80
  81
  82   // Reference tracks
  83
  84   refFile = new TFile(refFileName, "READ");
  85   TTree *treeRef = (TTree*)refFile->Get("TreeTR0_Sorted");
  86   TBranch *trkRef = treeRef->GetBranch( refTreeName[nRefPlane-1] );
  87   TClonesArray *trkRefs = new TClonesArray("AliTrackReference", 100);
  88   trkRef->SetAddress(&trkRefs);
  89
  90   Int_t ntracks = trkRef->GetEntries();
  91   cout << ntracks << endl;
  92
  93   // Tracks
  94
  95   TFile *trackFile = new TFile(trackFileName[nRefPlane-1]);
  96
  97   TTree *tracks = (TTree*)trackFile->Get(treeName[nRefPlane-1]);
  98   TBranch *recTracks = tracks->GetBranch("tracks");
  99   AliTPCtrack *track = 0;
 100   recTracks->SetAddress(&track);
 101
 102   for(Int_t t=0; t<tracks->GetEntries(); t++) {
 103
 104     cout << t << "\r";
 105
 106     recTracks->GetEvent(t);
 107
 108     Int_t lab = track->GetLabel();
 109     if (lab < 0 || lab > 100000) continue;
 110
 111     Int_t pdg = gAlice->Particle(lab)->GetPdgCode();
 112     treeRef->GetEvent(lab);
 113
 114     Int_t nEntrRef = trkRefs->GetEntries();
 115
 116     if (gAlice->Particle(lab)->P() < minP ) continue;
 117     if (gAlice->Particle(lab)->P() > maxP ) continue;
 118     if (gAlice->Particle(lab)->GetFirstMother() != -1) continue;
 119
 120     if (!(nEntrRef)) continue;
 121
 122     Double_t timeTrue, timeTrack;
 123     Double_t trueLength;
 124
 125     Int_t index =  nEntrRef-1;
 126     if (nRefPlane == 1) index = 0;
 127
 128     timeTrue = ((AliTrackReference*)(*trkRefs)[index])->GetTime() * 1e12;
 129     trueLength = ((AliTrackReference*)(*trkRefs)[index])->GetLength();
 130
 131     hLength->Fill(trueLength - track->GetIntegratedLength());
 132
 133     for(Int_t i=0; i<nTypes; i++) {
 134       timeTrack = track->GetIntegratedTime(codes[i]);
 135       if (abs(pdg) == codes[i]) hDiffSame[i]->Fill(timeTrue - timeTrack);
 136     }
 137   }
 138
 139   TCanvas *c = new TCanvas();
 140   c->SetCanvasSize(640, 800);
 141   c->SetWindowSize(650, 830);
 142   c->Divide(2,3);
 143
 144   c->cd(1);
 145   hLength->Draw();
 146   hLength->SetXTitle("True Length - Measured Length [cm]");
 147
 148   for(Int_t i=0; i<nTypes; i++) {
 149     c->cd(2+i);
 150     gPad->SetGridx();
 151     gPad->SetGridy();
 152     hDiffSame[i]->Draw();
 153     hDiffSame[i]->SetXTitle("time true - measured [ps]");
 154     hDiffSame[i]->Fit("gaus");
 155   }
 156
 157   char outFileName[80];
 158   const char *refPlaneNames[] = {"inTPC", "outTPC","outTRD"};
 159
 160   c->Modified();
 161   c->ForceUpdate();
 162
 163   sprintf(outFileName, "time_%s_%3.2f-%3.2f.ps", refPlaneNames[nRefPlane-1], minP, maxP);
 164   c->Print(outFileName);
 165
 166   if (treeRef) delete treeRef;
 167   if (tracks) delete tracks;
 168
 169 }