]>
Commit | Line | Data |
---|---|---|
b9de75e1 | 1 | /**************************************************************************** |
2 | * This macro is supposed to do reconstruction in the barrel ALICE trackers * | |
3 | * (Make sure you have TPC digits and ITS hits before using this macro !!!) * | |
4 | * It does the following steps (April 12, 2001): * | |
5 | * 1) TPC cluster finding * | |
6 | * 2) TPC track finding * | |
7 | * 3) ITS cluster finding V2 (via fast points !) * | |
8 | * 4) ITS track finding V2 * | |
9 | * 5) ITS back track propagation V2 * | |
10 | * 6) TPC back track propagation * | |
11 | * (Origin: I.Belikov, CERN, Jouri.Belikov@cern.ch) * | |
12 | ****************************************************************************/ | |
13 | ||
14 | #ifndef __CINT__ | |
15 | #include "alles.h" | |
16 | #include "AliMagF.h" | |
17 | #include "AliTPCtracker.h" | |
18 | ||
19 | #include "AliITS.h" | |
20 | #include "AliITSgeom.h" | |
21 | #include "AliITSRecPoint.h" | |
22 | #include "AliITSclusterV2.h" | |
23 | #include "AliITSsimulationFastPoints.h" | |
24 | #include "AliITStrackerV2.h" | |
25 | #endif | |
26 | ||
743a19f2 | 27 | Int_t TPCFindClusters(const Char_t *inname, const Char_t *outname, Int_t n); |
28 | Int_t TPCFindTracks(const Char_t *inname, const Char_t *outname, Int_t n); | |
29 | Int_t TPCSortTracks(const Char_t *inname, const Char_t *inname2, const Char_t *outname, Int_t n); | |
b9de75e1 | 30 | Int_t TPCPropagateBack(const Char_t *inname, const Char_t *outname); |
31 | ||
743a19f2 | 32 | Int_t ITSFindClusters(const Char_t *inname, const Char_t *outname, Int_t n); |
33 | Int_t ITSFindTracks(const Char_t *inname, const Char_t *inname2, const Char_t *outname, Int_t n); | |
b9de75e1 | 34 | Int_t ITSPropagateBack(const Char_t *inname, const Char_t *outname); |
35 | ||
36 | ||
6dcef50e | 37 | Int_t AliBarrelReconstruction(Int_t n=1) { |
106ea0fc | 38 | const Char_t *TPCdigName="rfio:galice.root"; |
b9de75e1 | 39 | const Char_t *TPCclsName="AliTPCclusters.root"; |
40 | const Char_t *TPCtrkName="AliTPCtracks.root"; | |
743a19f2 | 41 | const Char_t *TPCtrkNameS="AliTPCtracksSorted.root"; |
42 | ||
b9de75e1 | 43 | |
106ea0fc | 44 | const Char_t *ITSdigName="rfio:galice.root"; |
b9de75e1 | 45 | const Char_t *ITSclsName="AliITSclustersV2.root"; |
46 | const Char_t *ITStrkName="AliITStracksV2.root"; | |
47 | ||
48 | AliKalmanTrack::SetConvConst(100/0.299792458/0.2/gAlice->Field()->Factor()); | |
49 | ||
743a19f2 | 50 | |
b9de75e1 | 51 | // ********** Find TPC clusters *********** // |
743a19f2 | 52 | if (TPCFindClusters(TPCdigName,TPCclsName, n)) { |
b9de75e1 | 53 | cerr<<"Failed to get TPC clusters !\n"; |
54 | return 1; | |
55 | } | |
56 | ||
57 | // ********** Find TPC tracks *********** // | |
743a19f2 | 58 | if (TPCFindTracks(TPCclsName,TPCtrkName,n)) { |
b9de75e1 | 59 | cerr<<"Failed to get TPC tracks !\n"; |
60 | return 1; | |
743a19f2 | 61 | } |
62 | ||
b9de75e1 | 63 | // ********** Sort and label TPC tracks *********** // |
743a19f2 | 64 | if (TPCSortTracks(TPCclsName,TPCtrkName,TPCtrkNameS,n)) { |
b9de75e1 | 65 | cerr<<"Failed to sort TPC tracks !\n"; |
66 | return 1; | |
743a19f2 | 67 | } |
68 | ||
b9de75e1 | 69 | // ********** Find ITS clusters *********** // |
743a19f2 | 70 | if (ITSFindClusters(ITSdigName,ITSclsName,n)) { |
b9de75e1 | 71 | cerr<<"Failed to get ITS clusters !\n"; |
72 | return 1; | |
73 | } | |
743a19f2 | 74 | |
b9de75e1 | 75 | // ********** Find ITS tracks *********** // |
743a19f2 | 76 | { |
77 | //TFile *clsFile=TFile::Open(ITSclsName); | |
78 | if (ITSFindTracks(TPCtrkNameS,ITSclsName,ITStrkName,n)) { | |
b9de75e1 | 79 | cerr<<"Failed to get ITS tracks !\n"; |
80 | return 1; | |
81 | } | |
743a19f2 | 82 | //clsFile->Close(); |
83 | } | |
84 | return 1; | |
b9de75e1 | 85 | // ********** Back propagation of the ITS tracks *********** // |
86 | {TFile *clsFile=TFile::Open(ITSclsName); | |
743a19f2 | 87 | if (ITSPropagateBack(ITStrkName,TPCtrkNameS)) { |
b9de75e1 | 88 | cerr<<"Failed to propagate back the ITS tracks !\n"; |
89 | return 1; | |
90 | } | |
91 | clsFile->Close();} | |
92 | ||
93 | ||
94 | // ********** Back propagation of the TPC tracks *********** // | |
95 | {TFile *clsFile=TFile::Open(TPCclsName); | |
96 | if (TPCPropagateBack(TPCtrkName,TPCtrkName)) { | |
97 | cerr<<"Failed to propagate back the TPC tracks !\n"; | |
98 | return 1; | |
99 | } | |
100 | clsFile->Close();} | |
101 | ||
102 | return 0; | |
103 | } | |
104 | ||
105 | ||
743a19f2 | 106 | Int_t TPCFindClusters(const Char_t *inname, const Char_t *outname, Int_t n) { |
b9de75e1 | 107 | Int_t rc=0; |
108 | const Char_t *name="TPCFindClusters"; | |
109 | cerr<<'\n'<<name<<"...\n"; | |
110 | gBenchmark->Start(name); | |
111 | ||
112 | TFile *out=TFile::Open(outname,"recreate"); | |
113 | TFile *in =TFile::Open(inname); | |
114 | ||
115 | AliTPCParam *param=(AliTPCParam *)in->Get("75x40_100x60"); | |
116 | if (!param) {cerr<<"TPC parameters have not been found !\n"; return 1;} | |
117 | AliTPCv2 tpc; | |
118 | tpc.SetParam(param); | |
119 | ||
743a19f2 | 120 | //tpc.Digits2Clusters(out); //MI change |
121 | for (Int_t i=0;i<n;i++){ | |
122 | printf("Processing event %d\n",i); | |
123 | tpc.Digits2Clusters(out,i); | |
124 | // AliTPCclusterer::Digits2Clusters(dig, out, i); | |
125 | } | |
b9de75e1 | 126 | in->Close(); |
127 | out->Close(); | |
128 | gBenchmark->Stop(name); | |
129 | gBenchmark->Show(name); | |
130 | ||
131 | return rc; | |
132 | } | |
133 | ||
743a19f2 | 134 | Int_t TPCFindTracks(const Char_t *inname, const Char_t *outname, Int_t n) { |
b9de75e1 | 135 | Int_t rc=0; |
136 | const Char_t *name="TPCFindTracks"; | |
137 | cerr<<'\n'<<name<<"...\n"; | |
138 | gBenchmark->Start(name); | |
139 | TFile *out=TFile::Open(outname,"recreate"); | |
140 | TFile *in =TFile::Open(inname); | |
141 | AliTPCParam *param=(AliTPCParam *)in->Get("75x40_100x60"); | |
142 | if (!param) {cerr<<"TPC parameters have not been found !\n"; return 1;} | |
143 | ||
743a19f2 | 144 | //AliTPCtracker *tracker=new AliTPCtracker(param); |
145 | //rc=tracker->Clusters2Tracks(0,out); | |
146 | //delete tracker; | |
147 | ||
148 | for (Int_t i=0;i<n;i++){ | |
149 | printf("Processing event %d\n",i); | |
150 | AliTPCtracker *tracker = new AliTPCtracker(param,i); | |
151 | //Int_t rc= | |
152 | tracker->Clusters2Tracks(0,out); | |
153 | delete tracker; | |
154 | } | |
155 | ||
b9de75e1 | 156 | in->Close(); |
157 | out->Close(); | |
158 | gBenchmark->Stop(name); | |
159 | gBenchmark->Show(name); | |
160 | ||
161 | return rc; | |
162 | } | |
743a19f2 | 163 | Int_t TPCSortTracks(const Char_t *inname, const Char_t *inname2, const Char_t * outname, Int_t eventn){ |
b9de75e1 | 164 | Int_t rc=0; |
165 | const Char_t *name="TPCSortTracks"; | |
166 | cerr<<'\n'<<name<<"...\n"; | |
167 | gBenchmark->Start(name); | |
168 | ||
743a19f2 | 169 | TFile *out =TFile::Open(outname,"recreate"); |
170 | TFile *in1 =TFile::Open(inname); | |
171 | TFile *in2 =TFile::Open(inname2); | |
172 | ||
173 | ||
174 | AliTPCParam *param=(AliTPCParam *)in1->Get("75x40_100x60"); | |
b9de75e1 | 175 | if (!param) {cerr<<"TPC parameters have not been found !\n"; return 1;} |
176 | ||
b9de75e1 | 177 | |
743a19f2 | 178 | |
179 | // loop over events | |
180 | for (Int_t event=0;event<eventn; event++){ | |
181 | TObjArray tarray(10000); | |
182 | AliTPCtrack *iotrack=0; | |
183 | Int_t i; | |
184 | ||
185 | ||
186 | in2->cd(); | |
187 | char tname[100]; | |
188 | sprintf(tname,"TreeT_TPC_%d",event); | |
189 | ||
190 | TTree *tracktree=(TTree*)in2->Get(tname); | |
191 | TBranch *tbranch=tracktree->GetBranch("tracks"); | |
192 | Int_t nentr=(Int_t)tracktree->GetEntries(); | |
193 | for (i=0; i<nentr; i++) { | |
b9de75e1 | 194 | iotrack=new AliTPCtrack; |
195 | tbranch->SetAddress(&iotrack); | |
196 | tracktree->GetEvent(i); | |
197 | tarray.AddLast(iotrack); | |
743a19f2 | 198 | } |
199 | tarray.Sort(); | |
200 | // in2->Close(); | |
201 | ||
202 | //assign thacks GEANT labels | |
203 | in1->cd(); | |
204 | AliTPCtracker *tracker = new AliTPCtracker(param,event); | |
205 | tracker->LoadInnerSectors(); | |
206 | tracker->LoadOuterSectors(); | |
207 | for (i=0; i<nentr; i++) { | |
b9de75e1 | 208 | iotrack=(AliTPCtrack*)tarray.UncheckedAt(i); |
209 | tracker->CookLabel(iotrack,0.1); | |
743a19f2 | 210 | } |
211 | delete tracker; | |
212 | //in->Close(); | |
213 | //end of GEANT label assignment | |
214 | ||
215 | tracktree=new TTree(tname,"Tree with TPC tracks"); | |
216 | tracktree->Branch("tracks","AliTPCtrack",&iotrack,32000,0); | |
217 | for (i=0; i<nentr; i++) { | |
b9de75e1 | 218 | iotrack=(AliTPCtrack*)tarray.UncheckedAt(i); |
219 | tracktree->Fill(); | |
743a19f2 | 220 | } |
221 | out->cd(); | |
222 | tracktree->Write(); | |
b9de75e1 | 223 | } |
743a19f2 | 224 | |
b9de75e1 | 225 | out->Close(); |
743a19f2 | 226 | in2->Close(); |
227 | in1->Close(); | |
b9de75e1 | 228 | |
229 | gBenchmark->Stop(name); | |
230 | gBenchmark->Show(name); | |
231 | ||
232 | return rc; | |
233 | } | |
234 | ||
743a19f2 | 235 | Int_t ITSFindClusters(const Char_t *inname, const Char_t *outname, Int_t n) { |
b9de75e1 | 236 | Int_t rc=0; |
237 | const Char_t *name="ITSFindClusters"; | |
238 | cerr<<'\n'<<name<<"...\n"; | |
239 | gBenchmark->Start(name); | |
240 | TFile *out=TFile::Open(outname,"recreate"); | |
241 | TFile *in =TFile::Open(inname,"update"); | |
242 | ||
243 | if (!(gAlice=(AliRun*)in->Get("gAlice"))) { | |
244 | cerr<<"Can't get gAlice !\n"; | |
245 | return 1; | |
246 | } | |
743a19f2 | 247 | |
b9de75e1 | 248 | AliITS *ITS = (AliITS*)gAlice->GetModule("ITS"); |
249 | if (!ITS) { cerr<<"Can't get the ITS !\n"; return 1;} | |
b9de75e1 | 250 | AliITSgeom *geom=ITS->GetITSgeom(); |
743a19f2 | 251 | out->cd(); |
b9de75e1 | 252 | geom->Write(); |
743a19f2 | 253 | |
254 | Int_t ev=0; | |
255 | for (ev = 0; ev<n; ev++){ | |
256 | in->cd(); // !!!! MI directory must point to galice. - othervise problem with Tree -connection | |
257 | gAlice->GetEvent(ev); | |
258 | //gAlice->TreeR()->Reset(); //reset reconstructed tree | |
259 | ||
260 | ||
261 | TTree *pTree=gAlice->TreeR(); | |
262 | if (!pTree){ | |
263 | gAlice->MakeTree("R"); | |
264 | pTree = gAlice->TreeR(); | |
265 | } | |
266 | TBranch *branch=pTree->GetBranch("ITSRecPoints"); | |
267 | if (!branch) { | |
268 | //if not reconstructed ITS branch do reconstruction | |
269 | ITS->MakeBranch("R",0); | |
270 | //////////////// Taken from ITSHitsToFastPoints.C /////////////////////// | |
8d9c76b5 | 271 | for (Int_t i=0;i<3;i++) { |
272 | ITS->SetSimulationModel(i,new AliITSsimulationFastPoints()); | |
273 | } | |
743a19f2 | 274 | Int_t nsignal=25; |
275 | Int_t size=-1; | |
276 | Int_t bgr_ev=Int_t(ev/nsignal); | |
277 | ITS->HitsToFastRecPoints(ev,bgr_ev,size," ","All"," "); | |
278 | ////////////////////////////////////////////////////////////////////////// | |
279 | gAlice->GetEvent(ev); //MI comment - in HitsToFast... they reset treeR to 0 | |
280 | //they overwrite full reconstructed event ???? ... so lets connect TreeR one more | |
281 | //time | |
282 | } | |
283 | ||
284 | ||
285 | ||
286 | out->cd(); | |
287 | TClonesArray *clusters=new TClonesArray("AliITSclusterV2",10000); | |
288 | char cname[100]; | |
289 | sprintf(cname,"TreeC_ITS_%d",ev); | |
290 | ||
291 | TTree *cTree=new TTree(cname,"ITS clusters"); | |
292 | cTree->Branch("Clusters",&clusters); | |
293 | ||
294 | pTree=gAlice->TreeR(); | |
295 | if (!pTree) { cerr<<"Can't get TreeR !\n"; return 1; } | |
296 | branch=pTree->GetBranch("ITSRecPoints"); | |
297 | if (!branch) { cerr<<"Can't get ITSRecPoints branch !\n"; return 1;} | |
298 | TClonesArray *points=new TClonesArray("AliITSRecPoint",10000); | |
299 | branch->SetAddress(&points); | |
300 | ||
301 | TClonesArray &cl=*clusters; | |
302 | Int_t nclusters=0; | |
303 | Int_t nentr=(Int_t)pTree->GetEntries(); | |
304 | AliITSgeom *geom=ITS->GetITSgeom(); | |
305 | ||
306 | for (Int_t i=0; i<nentr; i++) { | |
b9de75e1 | 307 | if (!pTree->GetEvent(i)) {cTree->Fill(); continue;} |
308 | Int_t lay,lad,det; geom->GetModuleId(i,lay,lad,det); | |
309 | Float_t x,y,zshift; geom->GetTrans(lay,lad,det,x,y,zshift); | |
310 | Double_t rot[9]; geom->GetRotMatrix(lay,lad,det,rot); | |
311 | Double_t yshift = x*rot[0] + y*rot[1]; | |
312 | Int_t ndet=(lad-1)*geom->GetNdetectors(lay) + (det-1); | |
313 | Int_t ncl=points->GetEntriesFast(); | |
314 | nclusters+=ncl; | |
315 | for (Int_t j=0; j<ncl; j++) { | |
743a19f2 | 316 | AliITSRecPoint *p=(AliITSRecPoint*)points->UncheckedAt(j); |
317 | Float_t lp[5]; | |
318 | lp[0]=-p->GetX()-yshift; if (lay==1) lp[0]=-lp[0]; | |
319 | lp[1]=p->GetZ()+zshift; | |
320 | lp[2]=p->GetSigmaX2(); | |
321 | lp[3]=p->GetSigmaZ2(); | |
322 | lp[4]=p->GetQ(); | |
323 | Int_t lab[6]; | |
324 | lab[0]=p->GetLabel(0);lab[1]=p->GetLabel(1);lab[2]=p->GetLabel(2); | |
325 | lab[3]=ndet; | |
326 | ||
327 | Int_t label=lab[0]; | |
328 | TParticle *part=(TParticle*)gAlice->Particle(label); | |
329 | label=-3; | |
330 | while (part->P() < 0.005) { | |
331 | Int_t m=part->GetFirstMother(); | |
332 | if (m<0) {cerr<<"Primary momentum: "<<part->P()<<endl; break;} | |
333 | label=m; | |
334 | part=(TParticle*)gAlice->Particle(label); | |
335 | } | |
336 | if (lab[1]<0) lab[1]=label; | |
337 | else if (lab[2]<0) lab[2]=label; | |
338 | else cerr<<"No empty labels !\n"; | |
339 | ||
340 | new(cl[j]) AliITSclusterV2(lab,lp); | |
b9de75e1 | 341 | } |
342 | cTree->Fill(); clusters->Delete(); | |
343 | points->Delete(); | |
743a19f2 | 344 | } |
345 | cTree->Write(); | |
346 | cerr<<"Number of clusters: "<<nclusters<<endl; | |
347 | delete cTree; delete clusters; delete points; | |
b9de75e1 | 348 | |
743a19f2 | 349 | } |
b9de75e1 | 350 | |
351 | delete gAlice; gAlice=0; | |
352 | in->Close(); | |
353 | out->Close(); | |
354 | gBenchmark->Stop(name); | |
355 | gBenchmark->Show(name); | |
356 | ||
357 | return rc; | |
358 | } | |
359 | ||
743a19f2 | 360 | Int_t ITSFindTracks(const Char_t * inname, const Char_t *inname2, const Char_t *outname, Int_t n) { |
b9de75e1 | 361 | Int_t rc=0; |
362 | const Char_t *name="ITSFindTracks"; | |
363 | cerr<<'\n'<<name<<"...\n"; | |
364 | gBenchmark->Start(name); | |
365 | ||
b9de75e1 | 366 | |
367 | TFile *out=TFile::Open(outname,"recreate"); | |
368 | TFile *in =TFile::Open(inname); | |
743a19f2 | 369 | TFile *in2 =TFile::Open(inname2); |
370 | ||
371 | AliITSgeom *geom=(AliITSgeom*)gFile->Get("AliITSgeom"); | |
372 | if (!geom) { cerr<<"can't get ITS geometry !\n"; return 1;} | |
373 | ||
374 | ||
375 | for (Int_t i=0;i<n;i++){ | |
376 | AliITStrackerV2 tracker(geom,i); | |
377 | rc=tracker.Clusters2Tracks(in,out); | |
378 | } | |
379 | ||
b9de75e1 | 380 | in->Close(); |
743a19f2 | 381 | in2->Close(); |
b9de75e1 | 382 | out->Close(); |
383 | ||
384 | gBenchmark->Stop(name); | |
385 | gBenchmark->Show(name); | |
386 | ||
387 | return rc; | |
388 | } | |
389 | ||
390 | Int_t ITSPropagateBack(const Char_t *inname, const Char_t *outname) { | |
743a19f2 | 391 | |
392 | Int_t rc=0; | |
393 | /* | |
b9de75e1 | 394 | const Char_t *name="ITSPropagateBack"; |
395 | cerr<<'\n'<<name<<"...\n"; | |
396 | gBenchmark->Start(name); | |
397 | ||
398 | AliITSgeom *geom=(AliITSgeom*)gFile->Get("AliITSgeom"); | |
399 | if (!geom) { cerr<<"can't get ITS geometry !\n"; return 1;} | |
400 | AliITStrackerV2 tracker(geom); | |
401 | ||
402 | TFile *out=TFile::Open(outname,"update"); | |
403 | TFile *in =TFile::Open(inname); | |
404 | rc=tracker.PropagateBack(in,out); | |
405 | in->Close(); | |
406 | out->Close(); | |
407 | ||
408 | gBenchmark->Stop(name); | |
409 | gBenchmark->Show(name); | |
743a19f2 | 410 | */ |
b9de75e1 | 411 | return rc; |
412 | } | |
413 | ||
414 | Int_t TPCPropagateBack(const Char_t *inname, const Char_t *outname) { | |
415 | Int_t rc=0; | |
416 | const Char_t *name="TPCPropagateBack"; | |
417 | cerr<<'\n'<<name<<"...\n"; | |
418 | gBenchmark->Start(name); | |
419 | ||
420 | AliTPCParam *param=(AliTPCParam *)gFile->Get("75x40_100x60"); | |
421 | if (!param) {cerr<<"TPC parameters have not been found !\n"; return 1;} | |
422 | AliTPCtracker *tracker=new AliTPCtracker(param); | |
423 | ||
424 | TFile *out=TFile::Open(outname,"update"); | |
425 | TFile *in =TFile::Open(inname); | |
426 | rc=tracker->PropagateBack(in,out); | |
427 | delete tracker; | |
428 | in->Close(); | |
429 | out->Close(); | |
430 | ||
431 | gBenchmark->Stop(name); | |
432 | gBenchmark->Show(name); | |
433 | ||
434 | return rc; | |
435 | } | |
436 |