CMake: Adding missing ALI_DATE flag
[u/mrichter/AliRoot.git] / TPC / AliBarrelRec_TPCparam.C
CommitLineData
a7a1dd76 1/// \file AliBarrelRec_TPCparam.C
2/// \brief Macro to create AliESDs.root using parametrized TPC tracking and AliITStrackerSA (MI + SA)
3///
4/// Input files:
5/// - galice.root
6/// - Kinematics.root
7/// - TrackReferences.root
8/// - ITS.RecPoints.root (use AliRecontruction class)
9/// - ITS.Vertex.root (use $ALICE_ROOT/ITS/AliITSVertexerZTest.C)
10///
11/// \author A. Dainese - INFN Legnaro
056891c0 12
a7a1dd76 13void AliBarrelRec_TPCparam(Int_t firstEvent=0,Int_t lastEvent=0) {
64cf164e 14
970d08c3 15 Int_t collcode = 1; // pp collisions
225a9242 16 Bool_t useMeanVtx = kFALSE;
970d08c3 17
64cf164e 18 AliGeomManager::LoadGeometry("geometry.root");
970d08c3 19
20 if (gAlice) {
33c3c91a 21 delete AliRunLoader::Instance();
970d08c3 22 delete gAlice;
23 gAlice=0;
24 }
25 AliRunLoader *rl = AliRunLoader::Open("galice.root");
26 if (rl == 0x0) {
27 cerr<<"Can not open session"<<endl;
28 return;
29 }
30 Int_t retval = rl->LoadgAlice();
31 if (retval) {
32 cerr<<"LoadgAlice returned error"<<endl;
33 delete rl;
34 return;
35 }
36 retval = rl->LoadHeader();
37 if (retval) {
38 cerr<<"LoadHeader returned error"<<endl;
39 delete rl;
40 return;
41 }
42 gAlice=rl->GetAliRun();
43
44
45
46 // Get field from galice.root
47 AliMagF *fiel = (AliMagF*)gAlice->Field();
48 Double_t fieval=TMath::Abs((Double_t)fiel->SolenoidField()/10.);
49 // Set the conversion constant between curvature and Pt
50 AliTracker::SetFieldMap(fiel,kTRUE);
51
52 /**** The TPC corner ********************/
53
54 AliTPCtrackerParam tpcTrackerPar(collcode,fieval);
55 tpcTrackerPar.Init();
56
57 //**** Switch on the PID class (mandatory!)
58 AliPID pid;
59
60 /**** The ITS corner ********************/
61
62 AliITSLoader* itsl = (AliITSLoader*)rl->GetLoader("ITSLoader");
63 if (itsl == 0x0) {
64 cerr<<"Cannot get the ITS loader"<<endl;
65 return;
66 }
67 itsl->LoadRecPoints("read");
68
64cf164e 69 AliITSRecoParam * itsRecoParam = AliITSRecoParam::GetLowFluxParam();
70 AliITSReconstructor::SetRecoParam(itsRecoParam);
970d08c3 71
72 // Instance of the ITS tracker
64cf164e 73 AliITStrackerSA itsTracker(0);
970d08c3 74 Int_t ITSclusters[6] = {1,1,1,1,1,1};
64cf164e 75 itsTracker.SetLayersNotToSkip(ITSclusters);
970d08c3 76
77 // Primary vertex reconstruction in pp
78 AliESDVertex *initVertex = 0;
79 TFile *invtx = new TFile("AliESDVertexMean.root");
80 if(collcode==1 && useMeanVtx) {
81 // open the mean vertex
82 initVertex = (AliESDVertex*)invtx->Get("vtxmean");
83 } else {
84 Double_t pos[3]={0.5,0.5,0.};
85 Double_t err[3]={3.,3.,5.};
86 initVertex = new AliESDVertex(pos,err);
87 }
88 invtx->Close();
89 delete invtx;
64cf164e 90 AliVertexerTracks *vertexer = new AliVertexerTracks(AliTracker::GetBz());
970d08c3 91 vertexer->SetVtxStart(initVertex);
92 vertexer->SetDebug(0);
93 delete initVertex;
94 initVertex=0;
95
96 /***** The TREE for ESD is born *****/
97 TTree *esdTree=new TTree("esdTree","Tree with ESD objects");
64cf164e 98 AliESDEvent *event=0; AliESDEvent *eventTPCin=0;
99 event = new AliESDEvent();
100 event->CreateStdContent();
101 event->WriteToTree(esdTree);
970d08c3 102
103 if(firstEvent>rl->GetNumberOfEvents()) firstEvent=rl->GetNumberOfEvents()-1;
104 if(lastEvent>rl->GetNumberOfEvents()) lastEvent=rl->GetNumberOfEvents()-1;
105 cout<<" Number of events: "<<1+lastEvent-firstEvent<<endl;
106
107 TFile *ppZ = TFile::Open("ITS.Vertex.root"); // z vertices from SPD
108 AliESDVertex *vertexSPD = new AliESDVertex();
109 Char_t zver[100];
110 Double_t vtx[3]={0,0,0};
111 Double_t sigmavtx[3]={0.07,0.07,0.1};
112
113
114 //<---------------------------------- The Loop over events begins
115 TStopwatch timer;
116 Int_t trc;
117 for(Int_t i=firstEvent; i<=lastEvent; i++) {
118
225a9242 119 cout<<" Processing event number : "<<i<<endl;
64cf164e 120 //AliESDEvent *event = new AliESDEvent();
970d08c3 121 event->SetRunNumber(gAlice->GetRunNumber());
64cf164e 122 event->SetEventNumberInFile(i);
970d08c3 123 event->SetMagneticField(gAlice->Field()->SolenoidField());
124 rl->GetEvent(i);
125
126 //***** Primary vertex from SPD from file
127 sprintf(zver,"Event%d/Vertex",i);
128 vertexSPD = (AliESDVertex*)ppZ->Get(zver);
129 if(!vertexSPD) {
64cf164e 130 esdTree->Fill(); event->Reset();
970d08c3 131 continue;
132 }
133 event->SetVertex(vertexSPD);
134 vertexSPD->GetXYZ(vtx);
135 vertexSPD->GetSigmaXYZ(sigmavtx);
136
137 //***** TPC tracking
138 if ( (trc=tpcTrackerPar.BuildTPCtracks(event)) ) {
139 printf("exiting TPC tracker with code %d in event %d\n",trc,i);
64cf164e 140 esdTree->Fill(); event->Reset();
970d08c3 141 continue;
142 }
a654e257 143
225a9242 144 // make a copy of the ESD at this stage
145 eventTPCin = event;
146
970d08c3 147 //***** ITS tracking
148 itsTracker.AliTracker::SetVertex(vtx,sigmavtx);
64cf164e 149 // itsl->LoadRecPoints("read");
970d08c3 150 TTree *itsTree=itsl->TreeR();
151 if (!itsTree) {
152 cerr<<"Can't get the ITS cluster tree !\n";
64cf164e 153 esdTree->Fill(); event->Reset();
a654e257 154 return;
970d08c3 155 }
156 itsTracker.UnloadClusters();
157 itsTracker.LoadClusters(itsTree);
158 if ( (trc=itsTracker.Clusters2Tracks(event)) ) {
159 printf("exiting ITS tracker with code %d in event %d\n",trc,i);
64cf164e 160 esdTree->Fill(); event->Reset();
970d08c3 161 continue;
162 }
225a9242 163
164 // Bring kTPCin-tracks back to the TPC inner wall
165 BackToTPCInnerWall(event,eventTPCin);
166
167 // refit inward in ITS:
168 // - refit without vertex constraint
169 // - propagate through beam pipe to local x = 0
170 itsTracker.RefitInward(event);
970d08c3 171
970d08c3 172
173 //***** Vertex from ESD tracks
174 if(collcode==1) { // pp
175 AliESDVertex *vertexTrks =
176 (AliESDVertex*)vertexer->FindPrimaryVertex(event);
177 event->SetPrimaryVertex(vertexTrks);
178 }
179
180 esdTree->Fill();
64cf164e 181 event->Reset();
970d08c3 182
183 }//<-----------------------------------The Loop over events ends here
184 timer.Stop(); timer.Print();
185
186 // The AliESDs.root is born
187 TFile *ef = TFile::Open("AliESDs.root","RECREATE");
188 if (!ef || !ef->IsOpen()) {cerr<<"Can't open AliESDs.root !\n"; return;}
189
225a9242 190 //Write the tree and close everything
970d08c3 191 esdTree->Write();
192 delete esdTree;
193 ef->Close();
194
195 delete vertexer;
196 delete rl;
197
198 return;
056891c0 199}
4201c420 200//--------------------------------------------------------------------------
64cf164e 201void BackToTPCInnerWall(AliESDEvent *event,AliESDEvent *eventTPC) {
225a9242 202
203 Int_t ntracks = eventTPC->GetNumberOfTracks();
204 AliESDtrack *esdTrackTPC = 0;
4201c420 205
225a9242 206 // create relation between event and eventTPC
207 Int_t labelsTPC[100000000];
208 for(Int_t tr = 0; tr<ntracks; tr++) {
209 esdTrackTPC = (AliESDtrack*)event->GetTrack(tr);
210 labelsTPC[TMath::Abs(esdTrackTPC->GetLabel())] = tr;
211 }
212
213 ntracks = event->GetNumberOfTracks();
4201c420 214 AliESDtrack *esdTrack = 0;
225a9242 215 esdTrackTPC = 0;
216 Int_t indexTPC;
4201c420 217
218 // loop on tracks
225a9242 219 for(tr = 0; tr<ntracks; tr++) {
220 esdTrack = (AliESDtrack*)event->GetTrack(tr);
221 // set to kITSout the tracks that don't have kTPCin
222 // (they've been found by AliITStrackerSA)
223 if(!(esdTrack->GetStatus()&AliESDtrack::kTPCin)) {
224 esdTrack->SetStatus(AliESDtrack::kITSout);
225 continue;
226 }
227
228 // skip tracks that don't have kITSin
4201c420 229 if(!(esdTrack->GetStatus()&AliESDtrack::kITSin)) continue;
225a9242 230
231 indexTPC = labelsTPC[TMath::Abs(esdTrack->GetLabel())];
232 esdTrackTPC = (AliESDtrack*)eventTPC->GetTrack(indexTPC);
233
234 AliITStrackMI *itsTrack = 0;
4201c420 235 try {
225a9242 236 itsTrack = new AliITStrackMI(*esdTrackTPC);
4201c420 237 esdTrack = 0;
238 }
239 catch (const Char_t *msg) {
225a9242 240 Warning("ToTPCInnerWall",msg);
4201c420 241 continue;
242 }
225a9242 243 itsTrack->UpdateESDtrack(AliESDtrack::kITSout);
4201c420 244 esdTrack = new AliESDtrack(*(itsTrack->GetESDtrack()));
225a9242 245
4201c420 246 delete itsTrack;
247 } // end loop on tracks
248
970d08c3 249 return;
4201c420 250}
225a9242 251//--------------------------------------------------------------------------