New macros
[u/mrichter/AliRoot.git] / ITS / AliITSTracksV1.C
CommitLineData
af8e1c2d 1#ifndef __CINT__
2#include "AliITSgeom.h"
3#include "TParticle.h"
4 cerr<<"Reading tp good tracks...\n";
5 while (in>>gt[ngood].lab>>gt[ngood].code
6 >>gt[ngood].px >>gt[ngood].py>>gt[ngood].pz
7#endif
8
9struct GoodTrack {
10 Int_t fEventN;
11 Int_t lab;
12 Int_t code;
13 Bool_t flag;
14 Float_t px,py,pz;
15 Float_t x,y,z;
16 Float_t pxg,pyg,pzg,ptg;
17};
18
19
20Int_t AliITSTracksV1(Int_t evNumber1=0,Int_t evNumber2=0,Int_t nclust=5) {
21
22 cerr<<"Select tracks which have nclust rec points in ITS...\n";
23
24 GoodTrack gt[30000];
25 Int_t ngood=0;
26 ifstream in("good_tracks_tpc");
27 ofstream out("itsgood_tracks");
28
29 cerr<<"Reading good tracks...\n";
30 while (in>>gt[ngood].fEventN>>gt[ngood].lab>>gt[ngood].code
31 >>gt[ngood].px >>gt[ngood].py>>gt[ngood].pz
32 >>gt[ngood].x >>gt[ngood].y >>gt[ngood].z
33 >>gt[ngood].pxg >>gt[ngood].pyg >>gt[ngood].pzg
34 >>gt[ngood].ptg >>gt[ngood].flag) {
35 ngood++;
36 cerr<<ngood<<'\r';
37 if (ngood==30000) {
38 cerr<<"Too many good tracks !\n";
39 break;
40 }
41 }
42
43 if (!in.eof()) cerr<<"Read error (good_tracks_tpc) !\n";
44
45
46 TVector ntrackevtpc(evNumber2-evNumber1 +1);
47 Int_t nchange=-1;
48 Int_t nmev=-100;
49
50
51 Int_t i;
52 for( i =0; i<ngood ; i++){
53 if(gt[i].fEventN != nmev ){
54 nmev=gt[i].fEventN;
55 nchange++;
56 }
57 if(nmev == gt[i].fEventN) ntrackevtpc(nchange)++;
58 }
59
60
61 if (gClassTable->GetID("AliRun") < 0) {
62 gROOT->LoadMacro("loadlibs.C");
63 loadlibs();
64 } else {
65 delete gAlice;
66 gAlice=0;
67 }
68
69// Connect the Root Galice file containing Geometry, Kine and Hits
70 TFile *file = (TFile*)gROOT->GetListOfFiles()->FindObject("galice.root");
71 if (!file) file = new TFile("galice.root","UPDATE");
72 //if (!file) file = new TFile(filename);
73
74// Get AliRun object from file or create it if not on file
75 if (!gAlice) {
76 gAlice = (AliRun*)file->Get("gAlice");
77 if (gAlice) printf("AliRun object found on file\n");
78 if (!gAlice) gAlice = new AliRun("gAlice","Alice test program");
79 }
80
81
82 AliITS* ITS =(AliITS *)gAlice->GetDetector("ITS");
83 if (!ITS) return;
84
85 //TClonesArray *recPoints = ITS->RecPoints();
86 //cout<<" recPoints = "<<recPoints<<"\n";
87 Int_t numbpoints=0;
88 AliITSRecPoint *recp=0;
89 Int_t track=-3;
90 // Int_t modmin[]={0,80,240,324,500,1282};
91 // Int_t modmax[]={79,239,323,499,1281,2269};
92
93
94 Int_t modmin[6], modmax[6];
95 Int_t Nladder[6], Ndetector[6];
96 AliITSgeom *g1 = ITS->GetITSgeom();
97
98 for(Int_t i=0; i<6; i++) {
99 Nladder[i]=g1->GetNladders(i+1);
100 Ndetector[i]=g1->GetNdetectors(i+1);
101 //cout<<Nladder[i]<<" "<<Ndetector[i]<<"\n";
102 }
103 for(Int_t i=0; i<6; i++) {
104 modmin[i]=g1->GetModuleIndex(i+1,1,1);
105 modmax[i]=g1->GetModuleIndex(i+1,Nladder[i],Ndetector[i]);
106 //cout<<modmin[i]<<" "<<modmax[i]<<"\n";
107 }
108
109
110
111 for (int nev=evNumber1; nev<= evNumber2; nev++) {
112 Int_t nparticles = gAlice->GetEvent(nev);
113 cout << "nev " << nev <<endl;
114 cout << "nparticles " << nparticles <<endl;
115 if (nev < evNumber1) continue;
116 if (nparticles <= 0) return;
117 TClonesArray *recPoints = ITS->RecPoints();
118 TObjArray *particles=gAlice->Particles();
119 // TClonesArray *particles=gAlice->Particles();
120 // Int_t np=particles->GetEntriesFast();
121 Int_t np=gAlice->GetNtrack(); //FCA correction
122 TMatrix goodITS(np,6);
123 //Int_t modmin[]={0,80,240,324,500,1282};
124 //Int_t modmax[]={79,239,323,499,1281,2269};
125 // for(Int_t j=0; j<np; j++){
126 // for(Int_t j1=0; j1<6; j1++) goodITS(j,j1)=0;
127 // }
128
129
130 Int_t nent=gAlice->TreeR()->GetEntries();
131 printf("Found %d entries in the TreeR (must be one per module per event!)\n",nent); //
132
133 for (Int_t layer=1;layer<=6;layer++) {
134 for (Int_t mod=modmin[layer-1];mod<=modmax[layer-1];mod++) {
135 ITS->ResetRecPoints();
136 //gAlice->TreeR()->GetEvent(mod+1); //first entry in TreeR is empty
137 gAlice->TreeR()->GetEvent(mod); //modificato il 2-3-2001
138 numbpoints = recPoints->GetEntries();
139 if(!numbpoints) continue;
140 for (Int_t irec=0;irec<numbpoints;irec++) {
141 recp=(AliITSRecPoint*)recPoints->UncheckedAt(irec);
142 for (Int_t it=0;it<3;it++) {
143 track=recp->GetLabel(it);
144 if(track<0) continue;
145 if(track > np) {
146 cout<<" Error on track number \n";
147 getchar();
148 }
149 goodITS(track,layer-1)=1;
150 } //loop over tracks
151 } //loop over points
152 } //loop over modules
153 } //loop over layers
154
155 Int_t mingood=0, maxgood=0,jj=0;
156 if(nev==evNumber1) mingood=0;
157 else
158 {for(jj=evNumber1; jj<nev; jj++) mingood += ntrackevtpc(jj);}
159 for(jj=evNumber1; jj<=nev; jj++) maxgood+= ntrackevtpc(jj);
160 //cout<<" mingood maxgood = "<<mingood<<" "<<maxgood<<"\n";
161 for(i=mingood; i<maxgood; i++){
162 if(gt[i].lab>np) {cout<<" Error on TPC good tracks label \n"; getchar();}
163 Int_t ilab=gt[i].lab;
164 Int_t totclust=0;
165 for(Int_t j=0; j<6; j++) totclust+=goodITS(ilab,j);
166 if(totclust>=nclust && nev==gt[i].fEventN) {
167 gt[i].flag=1;
168 // printf("label nclusters %d %d\n",gt[i].lab,totclust); //
169 }
170 }
171
172 Int_t itsngood=0;
173
174 if (out) {
175 // for (Int_t ngd=0; ngd<ngood; ngd++) {
176 for(Int_t ngd=mingood; ngd<maxgood; ngd++){
177 TParticle *p = (TParticle*)gAlice->Particle(TMath::Abs(gt[ngd].lab)); //aggiunto 27-9
178 gt[ngd].x=p->Vx(); gt[ngd].y=p->Vy(); gt[ngd].z=p->Vz();
179 if(gt[ngd].flag) {
180 out<<gt[ngd].fEventN<<' '<<gt[ngd].lab<<' '<<gt[ngd].code<<' '
181 <<gt[ngd].px <<' '<<gt[ngd].py<<' '<<gt[ngd].pz<<' '
182 <<gt[ngd].x <<' '<<gt[ngd].y <<' '<<gt[ngd].z <<' '
183 <<gt[ngd].pxg <<' '<<gt[ngd].pyg <<' '<<gt[ngd].pzg <<' '
184 <<gt[ngd].ptg<<gt[ngd].flag<<endl;
185 itsngood++;
186 }
187 }
188 } else cerr<<"Can not open file (itsgood_tracks) !\n";
189
190
191 //cerr<<"Number of good tracks in TPC: "<<ngood<<endl;
192 cerr<<"Number of good tracks in TPC: "<<ntrackevtpc(nev)<<endl;
193 cerr<<"Number of good tracks in ITS: "<<itsngood<<endl;
194
195 } // event loop
196
197 out.close();
198 file->Close();
199
200 printf("after file close\n");
201
202 return 0;
203
204}
205