next50 trigger mask in AliHLTGlobalEsdConverterComponent
[u/mrichter/AliRoot.git] / TRD / Macros / AliTRDtoTOFanalysis.C
CommitLineData
59dfcadd 1#ifndef __CINT__
2 #include <iostream.h>
3
4 #include "AliTOF.h"
5 #include "AliTOFhit.h"
6 #include "AliTPCtrack.h"
7 #include "AliTRDtrack.h"
8 #include "AliTRDgeometry.h"
9
10 #include "alles.h"
11 #include "TFile.h"
12 #include "TParticle.h"
13 #include "TStopwatch.h"
14#endif
15
16void AliTRDtoTOFanalysis()
17{
18
19 Int_t nEvent = 0;
20 const Int_t maxIndex = 80000; // max number of primaries to be analysed
21 Int_t rtIndex[maxIndex];
22 for(Int_t i = 0; i < maxIndex; i++) rtIndex[i] = -1;
23
24 //****** tracking: Load reconstructed tracks
25
26 TFile *tf=TFile::Open("AliTRDtracks.root");
27
28 if (!tf->IsOpen()) {cerr<<"Can't open AliTRDtracks.root !\n"; return;}
29 TObjArray tarray(2000);
30
31 char tname[100];
32 sprintf(tname,"seedsTRDtoTOF2_%d",nEvent);
33 // for other radial positions use different trees:
34 // sprintf(tname,"seedsTRDtoTOF1_%d",nEvent);
35 // sprintf(tname,"seedsTRDtoPHOS_%d",nEvent);
36 // sprintf(tname,"seedsTRDtoRHIC_%d",nEvent);
37
38 TTree *tracktree=(TTree*)tf->Get(tname);
39
40 TBranch *tbranch=tracktree->GetBranch("tracks");
41
42 Int_t nRecTracks = (Int_t) tracktree->GetEntries();
43 cout<<"Found "<<nRecTracks<<" entries in the track tree "<<tname<<endl;
44
45 for (Int_t i=0; i<nRecTracks; i++) {
46 AliTPCtrack *iotrack=new AliTPCtrack();
47 tbranch->SetAddress(&iotrack);
48 tracktree->GetEvent(i);
49 tarray.AddLast(iotrack);
50 Int_t trackLabel = iotrack->GetLabel();
51
52 // printf("rt with %d clusters and label %d \n",
53 // iotrack->GetNumberOfClusters(), trackLabel);
54
55 if(trackLabel < 0) continue;
56 if(trackLabel >= maxIndex) continue;
57 rtIndex[trackLabel] = i;
58 }
59 tf->Close();
60
61 //********
62
63 TH1F *hdx = new TH1F("dx","X(hit) - X(track)",40,-10,10);
64 TH1F *hdy = new TH1F("dy","Y(hit) - Y(track)",100,-20,20);
65 TH1F *hdz = new TH1F("dz","Z(hit) - Z(track)",100,-20,20);
66 TH2F *hdt = new TH2F("dt","T(hit) vs S(track)",100,0,500,500,200,700);
67
68
69 // Connect the AliRoot file containing Geometry, Kine, Hits, and Digits
70 Char_t *alifile = "galice.root";
71
72 TFile *gafl = (TFile*) gROOT->GetListOfFiles()->FindObject(alifile);
73 if (!gafl) {
74 cout << "Open the ALIROOT-file " << alifile << endl;
75 gafl = new TFile(alifile);
76 }
77 else {
78 cout << alifile << " is already open" << endl;
79 }
80
81 // Get AliRun object from file or create it if not on file
82 gAlice = (AliRun*) gafl->Get("gAlice");
83 if (gAlice)
84 cout << "AliRun object found on file" << endl;
85 else
86 gAlice = new AliRun("gAlice","Alice test program");
87
88
89 // Import the Trees for the event nEvent in the file
90 const Int_t nparticles = gAlice->GetEvent(nEvent);
91 if (nparticles <= 0) return;
92
93 Float_t x,y,z,mass,e;
94 Double_t gX, gY, gZ; // rec. track coordinates in the global system
95 Double_t Px, Py, Pz; // rec. track momentum components in the global system
96 Float_t X, Y, Z; // local tracking coordinates for reconstructed track
97 Int_t nbytes = 0;
98 Int_t j,hit,ipart;
99 Int_t nhits;
100 Float_t tof;
101 TParticle *particle;
102 AliTPCtrack *rt;
103
104// Get pointers to Alice detectors and Hits containers
105 AliDetector *TOF = gAlice->GetDetector("TOF");
106 Int_t ntracks = (Int_t) gAlice->TreeH()->GetEntries();
107
108// Start loop on tracks in the hits containers
109
110 for (Int_t track=0; track < ntracks; track++) {
111
112 if(TOF) {
113 for(AliTOFhit* tofHit = (AliTOFhit*)TOF->FirstHit(track);
114 tofHit;
115 tofHit=(AliTOFhit*)TOF->NextHit()) {
116
117 ipart = tofHit->GetTrack();
118 if(ipart >= maxIndex) continue;
119 if(rtIndex[ipart] < 0) continue;
120
121 x = tofHit->X();
122 y = tofHit->Y();
123 z = tofHit->Z();
124
125 //******* tracking: extract track coordinates, momentum, etc.
126
127 rt = (AliTPCtrack*) tarray.UncheckedAt(rtIndex[ipart]);
128 Double_t tr_length = rt->GetLength(); // meaningfull only with ITS
129 // part fully implemented
130
131 AliTRDtrack *trd_rt = new AliTRDtrack(*rt, rt->GetAlpha());
132 trd_rt->GetGlobalXYZ(gX,gY,gZ); // returns running global coordinates
133 trd_rt->GetPxPyPz(Px,Py,Pz); // returns momentum in global system
134
135 /*
136 printf("\n hit position - rec. track position:\n");
137 printf(" X: %f - %f = %f \n", x, gX, x - gX);
138 printf(" Y: %f - %f = %f \n", y, gY, y - gY);
139 printf(" Z: %f - %f = %f \n", z, gZ, z - gZ);
140 printf("\n");
141 */
142
143 delete trd_rt;
144
145 // conversion from hit coordinates from global coordinate system
146 // to "local" tracking system
147
148 Float_t phi =TMath::ATan2(y,x);
149 if (phi > 2.*TMath::Pi()) phi -= 2.*TMath::Pi();
150 if (phi < 0. ) phi += 2.*TMath::Pi();
151 Int_t sec = Int_t(phi/AliTRDgeometry::GetAlpha()) % 18;
152 Float_t alpha = (sec + 0.5) * AliTRDgeometry::GetAlpha();
153
154 Float_t tmp=x*TMath::Cos(alpha) + y*TMath::Sin(alpha);
155 y= - x*TMath::Sin(alpha) + y*TMath::Cos(alpha);
156 x=tmp;
157
158 // particle = (TParticle*)gAlice->Particle(ipart);
159 // if (particle->GetFirstMother() < 0) continue;
160
161 X = (Float_t) rt->GetX(); // radial position in the tracking system
162 Y = (Float_t) rt->GetY(); // r*phi position within the sector
163 Z = (Float_t) rt->GetZ(); // Z position
164
165
166 if(TMath::Abs(X-x-1) > 2) continue;
167
168 hdx->Fill(X-x);
169 hdy->Fill(Y-y);
170 hdz->Fill(Z-z);
171
172 //**********
173
174 }
175 }
176 }
177
178 TCanvas* c1 = new TCanvas("c1", "c1", 210, 210, 910, 940);
179 c1->SetFillColor(10);
180 c1->Divide(2,2);
181 c1->cd(1); hdx->Draw();
182 c1->cd(2); hdy->Draw();
183 c1->cd(3); hdz->Draw();
184 c1->cd(4); hdt->Draw();
185}