]>
Commit | Line | Data |
---|---|---|
64bb243d | 1 | //******************************************************************** |
2 | // Example of the reconstruction that generates the ESD | |
3 | // Input files: | |
4 | // a) file containing the ITS clusters | |
5 | // (the AliITSFindClustersV2.C macro can be used to generate it) | |
6 | // b) file containing the TPC clusters | |
7 | // (the AliTPCFindClusters.C macro can be used to generate it) | |
8 | // c) file containing the TRD clusters | |
9 | // (the AliTRDdigits2cluster.C macro can be used to generate it) | |
10 | // d) file containing the TOF digits | |
11 | // (the AliTOFSDigits2Digits.C macro can be used to generate it) | |
12 | // Ouput file: | |
13 | // AliESDs.root containing the ESD events | |
14 | // | |
15 | // Origin: Iouri Belikov, CERN, Jouri.Belikov@cern.ch | |
16 | //******************************************************************** | |
17 | ||
18 | #if !defined(__CINT__) || defined(__MAKECINT__) | |
19 | #include <Riostream.h> | |
20 | #include "TFile.h" | |
21 | #include "TSystem.h" | |
22 | #include "TStopwatch.h" | |
23 | #include "TArrayF.h" | |
24 | ||
25 | #include "AliMagF.h" | |
26 | #include "AliRun.h" | |
27 | #include "AliRunLoader.h" | |
28 | #include "AliLoader.h" | |
29 | #include "AliHeader.h" | |
30 | #include "AliGenEventHeader.h" | |
31 | ||
32 | #include "AliESD.h" | |
33 | #include "AliESDpid.h" | |
34 | ||
35 | #include "AliITS.h" | |
36 | #include "AliITSgeom.h" | |
37 | #include "AliITStrackerV2.h" | |
38 | #include "AliV0vertexer.h" | |
39 | #include "AliCascadeVertexer.h" | |
40 | #include "AliITSpidESD.h" | |
41 | #include "AliITSLoader.h" | |
42 | ||
43 | #include "AliTPCParam.h" | |
44 | #include "AliTPCtracker.h" | |
45 | #include "AliTPCtrackerMI.h" | |
46 | #include "AliTPCpidESD.h" | |
47 | #include "AliTPCLoader.h" | |
48 | ||
49 | #include "AliTRDtracker.h" | |
50 | #include "AliTRDPartID.h" | |
51 | ||
52 | #include "AliTOFpidESD.h" | |
53 | #include "AliTOF.h" | |
54 | #include "AliTOFGeometry.h" | |
55 | #endif | |
56 | ||
57 | extern AliRun *gAlice; | |
58 | extern TFile *gFile; | |
59 | ||
60 | Int_t AliESDtest(Int_t nev=1,Int_t run=0, Int_t debug =1) { | |
61 | ||
62 | /**** Initialization of the NewIO *******/ | |
63 | ||
64 | if (gAlice) { | |
65 | delete gAlice->GetRunLoader(); | |
66 | delete gAlice; | |
67 | gAlice=0; | |
68 | } | |
69 | ||
70 | AliRunLoader *rl = AliRunLoader::Open("galice.root"); | |
71 | if (rl == 0x0) { | |
72 | cerr<<"Can not open session"<<endl; | |
73 | return 1; | |
74 | } | |
75 | Int_t retval = rl->LoadgAlice(); | |
76 | if (retval) { | |
77 | cerr<<"AliESDtest.C : LoadgAlice returned error"<<endl; | |
78 | delete rl; | |
79 | return 1; | |
80 | } | |
81 | retval = rl->LoadHeader(); | |
82 | if (retval) { | |
83 | cerr<<"AliESDtest.C : LoadHeader returned error"<<endl; | |
84 | delete rl; | |
85 | return 2; | |
86 | } | |
87 | gAlice=rl->GetAliRun(); | |
88 | ||
89 | ||
90 | AliKalmanTrack::SetConvConst( | |
91 | 1000/0.299792458/gAlice->Field()->SolenoidField() | |
92 | ); | |
93 | ||
94 | ||
95 | /**** The ITS corner ********************/ | |
96 | ||
97 | AliITSLoader* itsl = (AliITSLoader*)rl->GetLoader("ITSLoader"); | |
98 | if (itsl == 0x0) { | |
99 | cerr<<"AliESDtest.C : Can not get the ITS loader"<<endl; | |
100 | return 3; | |
101 | } | |
102 | itsl->LoadRecPoints("read"); | |
103 | ||
104 | AliITS *dITS = (AliITS*)gAlice->GetDetector("ITS"); | |
105 | if (!dITS) { | |
106 | cerr<<"AliESDtest.C : Can not find the ITS detector !"<<endl; | |
107 | return 4; | |
108 | } | |
109 | AliITSgeom *geom = dITS->GetITSgeom(); | |
110 | ||
111 | //An instance of the ITS tracker | |
112 | AliITStrackerV2 itsTracker(geom); | |
113 | ||
114 | //An instance of the ITS PID maker | |
115 | Double_t parITS[]={35.5,0.11,10.}; | |
116 | AliITSpidESD itsPID(parITS); | |
117 | ||
118 | //An instance of the V0 finder | |
119 | Double_t cuts[]={33, // max. allowed chi2 | |
120 | 0.16,// min. allowed negative daughter's impact parameter | |
121 | 0.05,// min. allowed positive daughter's impact parameter | |
122 | 0.080,// max. allowed DCA between the daughter tracks | |
123 | 0.998,// max. allowed cosine of V0's pointing angle | |
124 | 0.9, // min. radius of the fiducial volume | |
125 | 2.9 // max. radius of the fiducial volume | |
126 | }; | |
127 | AliV0vertexer vtxer(cuts); | |
128 | ||
129 | Double_t cts[]={33., // max. allowed chi2 | |
130 | 0.05, // min. allowed V0 impact parameter | |
131 | 0.008, // window around the Lambda mass | |
132 | 0.035, // min. allowed bachelor's impact parameter | |
133 | 0.10, // max. allowed DCA between a V0 and a track | |
134 | 0.9985, //max. allowed cosine of the cascade pointing angle | |
135 | 0.9, // min. radius of the fiducial volume | |
136 | 2.9 // max. radius of the fiducial volume | |
137 | }; | |
138 | AliCascadeVertexer cvtxer=AliCascadeVertexer(cts); | |
139 | ||
140 | /**** The TPC corner ********************/ | |
141 | AliDataLoader* dl1,*dl2,*dl3; | |
142 | AliTPCLoader* tpcl = (AliTPCLoader*)rl->GetLoader("TPCLoader"); | |
143 | if (tpcl == 0x0) { | |
144 | cerr<<"AliESDtest.C : can not get the TPC loader"<<endl; | |
145 | return 5; | |
146 | } | |
147 | tpcl->LoadRecPoints("read"); | |
148 | ||
149 | rl->CdGAFile(); | |
150 | AliTPCParam *par=(AliTPCParam *)gDirectory->Get("75x40_100x60_150x60"); | |
151 | if (!par) { | |
152 | cerr<<"TPC parameters have not been found !\n"; | |
153 | return 6; | |
154 | } | |
155 | // Debug option - making loaders tracks of TPC | |
156 | if (debug>0){ | |
157 | dl1 = new AliDataLoader("TPC.TracksOne.root", "TreeTOne", "TracksOne"); | |
158 | dl2 = new AliDataLoader("TPC.TracksTwo.root", "TreeTTwo", "TracksTwo"); | |
159 | dl3 = new AliDataLoader("TPC.TracksThree.root","TreeTThree", "TracksThree"); | |
160 | tpcl->AddDataLoader(dl1); | |
161 | tpcl->AddDataLoader(dl2); | |
162 | tpcl->AddDataLoader(dl3); | |
163 | } | |
164 | ||
165 | ||
166 | //An instance of the TPC tracker | |
167 | AliTPCtrackerMI tpcTracker(par); | |
168 | tpcTracker.SetDebug(debug); | |
169 | ||
170 | //An instance of the TPC PID maker | |
171 | Double_t parTPC[]={45.0,0.08,10.}; | |
172 | AliTPCpidESD tpcPID(parTPC); | |
173 | ||
174 | ||
175 | /**** The TRD corner ********************/ | |
176 | ||
177 | AliLoader* trdl = rl->GetLoader("TRDLoader"); | |
178 | if (trdl == 0x0) { | |
179 | cerr<<"AliESDtest.C : can not get the TRD loader"<<endl; | |
180 | return 5; | |
181 | } | |
182 | trdl->LoadRecPoints("read"); | |
183 | ||
184 | //An instance of the TRD tracker | |
185 | rl->CdGAFile(); | |
186 | AliTRDtracker trdTracker(gFile); //galice.root file | |
187 | ||
188 | /* | |
189 | //An instance of the TRD PID maker | |
190 | TFile* pidFile = TFile::Open("pid.root"); | |
191 | if (!pidFile->IsOpen()) { | |
192 | cerr << "Can't get pid.root !\n"; | |
193 | return 7; | |
194 | } | |
195 | AliTRDPartID* trdPID = (AliTRDPartID*) pidFile->Get("AliTRDPartID"); | |
196 | if (!trdPID) { | |
197 | cerr << "Can't get PID object !\n"; | |
198 | return 8; | |
199 | } | |
200 | */ | |
201 | ||
202 | ||
203 | /**** The TOF corner ********************/ | |
204 | AliTOF *dTOF = (AliTOF*)gAlice->GetDetector("TOF"); | |
205 | if (!dTOF) { | |
206 | cerr<<"AliESDtest.C : Can not find the TOF detector !"<<endl; | |
207 | return 4; | |
208 | } | |
209 | AliTOFGeometry *tofGeo = dTOF->GetGeometry(); | |
210 | if (!tofGeo) { | |
211 | cerr<<"AliESDtest.C : Can not find the TOF geometry !"<<endl; | |
212 | return 4; | |
213 | } | |
214 | ||
215 | AliLoader* tofl = rl->GetLoader("TOFLoader"); | |
216 | if (tofl == 0x0) { | |
217 | cerr<<"AliESDtest.C : can not get the TOF loader"<<endl; | |
218 | return 5; | |
219 | } | |
220 | tofl->LoadDigits("read"); | |
221 | ||
222 | //Instance of the TOF PID maker | |
223 | Double_t parTOF[]={130.,5.}; | |
224 | AliTOFpidESD tofPID(parTOF); | |
225 | ||
226 | ||
227 | //rl->UnloadgAlice(); | |
228 | ||
229 | ||
230 | TFile *ef=TFile::Open("AliESDs.root","RECREATE"); | |
231 | if (!ef->IsOpen()) {cerr<<"Can't AliESDs.root !\n"; return 1;} | |
232 | ||
233 | TStopwatch timer; | |
234 | Int_t rc=0; | |
235 | if (nev>rl->GetNumberOfEvents()) nev=rl->GetNumberOfEvents(); | |
236 | //The loop over events | |
237 | for (Int_t i=0; i<nev; i++) { | |
238 | cerr<<"\n\nProcessing event number : "<<i<<endl; | |
239 | AliESD *event=new AliESD(); | |
240 | event->SetRunNumber(run); | |
241 | event->SetEventNumber(i); | |
242 | ||
243 | rl->GetEvent(i); | |
244 | ||
245 | //***** Primary vertex reconstruction (MC vertex position, for the moment) | |
246 | TArrayF v(3); | |
247 | rl->GetHeader()->GenEventHeader()->PrimaryVertex(v); | |
248 | Double_t vtx[3]={v[0],v[1],v[2]}; | |
249 | Double_t cvtx[6]={ | |
250 | 0.005, | |
251 | 0.000, 0.005, | |
252 | 0.000, 0.000, 0.010 | |
253 | }; | |
254 | event->SetVertex(vtx,cvtx); | |
255 | cvtx[1]=cvtx[0]; cvtx[2]=cvtx[5]; //trackers use only the diag.elements | |
256 | ||
257 | //***** Initial path towards the primary vertex | |
258 | tpcTracker.SetVertex(vtx,cvtx); | |
259 | //printf("VERTEX- vtx,cvtx\n",vtx,cvtx); | |
260 | TTree *tpcTree=tpcl->TreeR(); | |
261 | if (!tpcTree) { | |
262 | cerr<<"Can't get the TPC cluster tree !\n"; | |
263 | return 4; | |
264 | } | |
265 | tpcTracker.LoadClusters(tpcTree); | |
266 | rc+=tpcTracker.Clusters2Tracks(event); | |
267 | if (debug){ | |
268 | dl1->Load("recreate"); | |
269 | tpcTracker.WriteTracks(dl1->Tree()); | |
270 | dl1->WriteData("RECREATE"); | |
271 | } | |
272 | tpcTracker.DeleteSeeds(); | |
273 | // | |
274 | ||
275 | tpcPID.MakePID(event); // preliminary PID | |
276 | AliESDpid::MakePID(event); // for the ITS tracker | |
277 | ||
278 | itsTracker.SetVertex(vtx,cvtx); | |
279 | TTree *itsTree=itsl->TreeR(); | |
280 | if (!itsTree) { | |
281 | cerr<<"Can't get the ITS cluster tree !\n"; | |
282 | return 4; | |
283 | } | |
284 | itsTracker.LoadClusters(itsTree); | |
285 | rc+=itsTracker.Clusters2Tracks(event); | |
286 | ||
287 | ||
288 | //***** Back propagation towards the outer barrel detectors | |
289 | rc+=itsTracker.PropagateBack(event); | |
290 | ||
291 | rc+=tpcTracker.PropagateBack(event); | |
292 | if (debug){ | |
293 | dl2->Load("recreate"); | |
294 | tpcTracker.WriteTracks(dl2->Tree()); | |
295 | dl2->WriteData("RECREATE"); | |
296 | } | |
297 | tpcTracker.DeleteSeeds(); | |
298 | ||
299 | TTree *trdTree=trdl->TreeR(); | |
300 | if (!trdTree) { | |
301 | cerr<<"Can't get the TRD cluster tree !\n"; | |
302 | return 4; | |
303 | } | |
304 | trdTracker.LoadClusters(trdTree); | |
305 | rc+=trdTracker.PropagateBack(event); | |
306 | ||
307 | ||
308 | /* | |
309 | for (Int_t iTrack = 0; iTrack < event->GetNumberOfTracks(); iTrack++) { | |
310 | AliESDtrack* track = event->GetTrack(iTrack); | |
311 | trdPID->MakePID(track); | |
312 | } | |
313 | */ | |
314 | ||
315 | TTree *tofTree=tofl->TreeD(); | |
316 | if (!tofTree) { | |
317 | cerr<<"Can't get the TOF cluster tree !\n"; | |
318 | return 4; | |
319 | } | |
320 | //tofPID.LoadClusters(tofTree,tofGeo); | |
321 | //tofPID.MakePID(event); | |
322 | //tofPID.UnloadClusters(); | |
323 | ||
324 | ||
325 | //***** Now the final refit at the primary vertex... | |
326 | rc+=trdTracker.RefitInward(event); | |
327 | trdTracker.UnloadClusters(); | |
328 | ||
329 | rc+=tpcTracker.RefitInward(event); | |
330 | if (debug){ | |
331 | dl3->Load("recreate"); | |
332 | tpcTracker.WriteTracks(dl3->Tree()); | |
333 | dl3->WriteData("RECREATE"); | |
334 | } | |
335 | tpcTracker.DeleteSeeds(); | |
336 | ||
337 | ||
338 | tpcTracker.UnloadClusters(); | |
339 | tpcPID.MakePID(event); | |
340 | ||
341 | rc+=itsTracker.RefitInward(event); | |
342 | itsTracker.UnloadClusters(); | |
343 | itsPID.MakePID(event); | |
344 | ||
345 | ||
346 | //***** Here is the combined PID | |
347 | AliESDpid::MakePID(event); | |
348 | ||
349 | ||
350 | //***** Hyperon reconstruction | |
351 | vtxer.SetVertex(vtx); | |
352 | rc+=vtxer.Tracks2V0vertices(event); // V0 finding | |
353 | rc+=cvtxer.V0sTracks2CascadeVertices(event); // cascade finding | |
354 | ||
355 | ||
356 | ||
357 | //***** Some final manipulations with this event | |
358 | if (rc==0) { | |
359 | Char_t ename[100]; | |
360 | sprintf(ename,"%d",i); | |
361 | ef->cd(); | |
362 | if (!event->Write(ename)) rc++; | |
363 | } | |
364 | if (rc) { | |
365 | cerr<<"Something bad happened...\n"; | |
366 | } | |
367 | delete event; | |
368 | } | |
369 | timer.Stop(); timer.Print(); | |
370 | ||
371 | //pidFile->Close(); | |
372 | ||
373 | delete par; | |
374 | ||
375 | ef->Close(); | |
376 | ||
377 | delete rl; | |
378 | ||
379 | return rc; | |
380 | } |