1 //********************************************************************
2 // Example of the reconstruction that generates the ESD
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)
13 // AliESDs.root containing the ESD events
15 // Origin: Iouri Belikov, CERN, Jouri.Belikov@cern.ch
16 //********************************************************************
18 #if !defined(__CINT__) || defined(__MAKECINT__)
19 #include <Riostream.h>
22 #include "TStopwatch.h"
27 #include "AliRunLoader.h"
28 #include "AliLoader.h"
29 #include "AliHeader.h"
30 #include "AliGenEventHeader.h"
33 #include "AliESDpid.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"
43 #include "AliTPCParam.h"
44 #include "AliTPCtracker.h"
45 #include "AliTPCtrackerMI.h"
46 #include "AliTPCpidESD.h"
47 #include "AliTPCLoader.h"
49 #include "AliTRDtracker.h"
50 #include "AliTRDPartID.h"
52 #include "AliTOFpidESD.h"
54 #include "AliTOFGeometry.h"
57 extern AliRun *gAlice;
60 Int_t AliESDtest(Int_t nev=1,Int_t run=0, Int_t debug =1) {
62 /**** Initialization of the NewIO *******/
65 delete gAlice->GetRunLoader();
70 AliRunLoader *rl = AliRunLoader::Open("galice.root");
72 cerr<<"Can not open session"<<endl;
75 Int_t retval = rl->LoadgAlice();
77 cerr<<"AliESDtest.C : LoadgAlice returned error"<<endl;
81 retval = rl->LoadHeader();
83 cerr<<"AliESDtest.C : LoadHeader returned error"<<endl;
87 gAlice=rl->GetAliRun();
90 AliKalmanTrack::SetConvConst(
91 1000/0.299792458/gAlice->Field()->SolenoidField()
95 /**** The ITS corner ********************/
97 AliITSLoader* itsl = (AliITSLoader*)rl->GetLoader("ITSLoader");
99 cerr<<"AliESDtest.C : Can not get the ITS loader"<<endl;
102 itsl->LoadRecPoints("read");
104 AliITS *dITS = (AliITS*)gAlice->GetDetector("ITS");
106 cerr<<"AliESDtest.C : Can not find the ITS detector !"<<endl;
109 AliITSgeom *geom = dITS->GetITSgeom();
111 //An instance of the ITS tracker
112 AliITStrackerV2 itsTracker(geom);
114 //An instance of the ITS PID maker
115 Double_t parITS[]={35.5,0.11,10.};
116 AliITSpidESD itsPID(parITS);
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
127 AliV0vertexer vtxer(cuts);
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
138 AliCascadeVertexer cvtxer=AliCascadeVertexer(cts);
140 /**** The TPC corner ********************/
141 AliDataLoader* dl1,*dl2,*dl3;
142 AliTPCLoader* tpcl = (AliTPCLoader*)rl->GetLoader("TPCLoader");
144 cerr<<"AliESDtest.C : can not get the TPC loader"<<endl;
147 tpcl->LoadRecPoints("read");
150 AliTPCParam *par=(AliTPCParam *)gDirectory->Get("75x40_100x60_150x60");
152 cerr<<"TPC parameters have not been found !\n";
155 // Debug option - making loaders tracks of TPC
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);
166 //An instance of the TPC tracker
167 AliTPCtrackerMI tpcTracker(par);
168 tpcTracker.SetDebug(debug);
170 //An instance of the TPC PID maker
171 Double_t parTPC[]={45.0,0.08,10.};
172 AliTPCpidESD tpcPID(parTPC);
175 /**** The TRD corner ********************/
177 AliLoader* trdl = rl->GetLoader("TRDLoader");
179 cerr<<"AliESDtest.C : can not get the TRD loader"<<endl;
182 trdl->LoadRecPoints("read");
184 //An instance of the TRD tracker
186 AliTRDtracker trdTracker(gFile); //galice.root file
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";
195 AliTRDPartID* trdPID = (AliTRDPartID*) pidFile->Get("AliTRDPartID");
197 cerr << "Can't get PID object !\n";
203 /**** The TOF corner ********************/
204 AliTOF *dTOF = (AliTOF*)gAlice->GetDetector("TOF");
206 cerr<<"AliESDtest.C : Can not find the TOF detector !"<<endl;
209 AliTOFGeometry *tofGeo = dTOF->GetGeometry();
211 cerr<<"AliESDtest.C : Can not find the TOF geometry !"<<endl;
215 AliLoader* tofl = rl->GetLoader("TOFLoader");
217 cerr<<"AliESDtest.C : can not get the TOF loader"<<endl;
220 tofl->LoadDigits("read");
222 //Instance of the TOF PID maker
223 Double_t parTOF[]={130.,5.};
224 AliTOFpidESD tofPID(parTOF);
227 //rl->UnloadgAlice();
230 TFile *ef=TFile::Open("AliESDs.root","RECREATE");
231 if (!ef->IsOpen()) {cerr<<"Can't AliESDs.root !\n"; return 1;}
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);
245 //***** Primary vertex reconstruction (MC vertex position, for the moment)
247 rl->GetHeader()->GenEventHeader()->PrimaryVertex(v);
248 Double_t vtx[3]={v[0],v[1],v[2]};
254 event->SetVertex(vtx,cvtx);
255 cvtx[1]=cvtx[0]; cvtx[2]=cvtx[5]; //trackers use only the diag.elements
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();
262 cerr<<"Can't get the TPC cluster tree !\n";
265 tpcTracker.LoadClusters(tpcTree);
266 rc+=tpcTracker.Clusters2Tracks(event);
268 dl1->Load("recreate");
269 tpcTracker.WriteTracks(dl1->Tree());
270 dl1->WriteData("RECREATE");
272 tpcTracker.DeleteSeeds();
275 tpcPID.MakePID(event); // preliminary PID
276 AliESDpid::MakePID(event); // for the ITS tracker
278 itsTracker.SetVertex(vtx,cvtx);
279 TTree *itsTree=itsl->TreeR();
281 cerr<<"Can't get the ITS cluster tree !\n";
284 itsTracker.LoadClusters(itsTree);
285 rc+=itsTracker.Clusters2Tracks(event);
288 //***** Back propagation towards the outer barrel detectors
289 rc+=itsTracker.PropagateBack(event);
291 rc+=tpcTracker.PropagateBack(event);
293 dl2->Load("recreate");
294 tpcTracker.WriteTracks(dl2->Tree());
295 dl2->WriteData("RECREATE");
297 tpcTracker.DeleteSeeds();
299 TTree *trdTree=trdl->TreeR();
301 cerr<<"Can't get the TRD cluster tree !\n";
304 trdTracker.LoadClusters(trdTree);
305 rc+=trdTracker.PropagateBack(event);
309 for (Int_t iTrack = 0; iTrack < event->GetNumberOfTracks(); iTrack++) {
310 AliESDtrack* track = event->GetTrack(iTrack);
311 trdPID->MakePID(track);
315 TTree *tofTree=tofl->TreeD();
317 cerr<<"Can't get the TOF cluster tree !\n";
320 //tofPID.LoadClusters(tofTree,tofGeo);
321 //tofPID.MakePID(event);
322 //tofPID.UnloadClusters();
325 //***** Now the final refit at the primary vertex...
326 rc+=trdTracker.RefitInward(event);
327 trdTracker.UnloadClusters();
329 rc+=tpcTracker.RefitInward(event);
331 dl3->Load("recreate");
332 tpcTracker.WriteTracks(dl3->Tree());
333 dl3->WriteData("RECREATE");
335 tpcTracker.DeleteSeeds();
338 tpcTracker.UnloadClusters();
339 tpcPID.MakePID(event);
341 rc+=itsTracker.RefitInward(event);
342 itsTracker.UnloadClusters();
343 itsPID.MakePID(event);
346 //***** Here is the combined PID
347 AliESDpid::MakePID(event);
350 //***** Hyperon reconstruction
351 vtxer.SetVertex(vtx);
352 rc+=vtxer.Tracks2V0vertices(event); // V0 finding
353 rc+=cvtxer.V0sTracks2CascadeVertices(event); // cascade finding
357 //***** Some final manipulations with this event
360 sprintf(ename,"%d",i);
362 if (!event->Write(ename)) rc++;
365 cerr<<"Something bad happened...\n";
369 timer.Stop(); timer.Print();