Changes to compile with Root6 on macosx64
[u/mrichter/AliRoot.git] / PWGPP / TPC / AliPerformanceTPC.cxx
1 //------------------------------------------------------------------------------
2 // Implementation of AliPerformanceTPC class. It keeps information from 
3 // comparison of reconstructed and MC particle tracks. In addtion, 
4 // it keeps selection cuts used during comparison. The comparison 
5 // information is stored in the ROOT histograms. Analysis of these 
6 // histograms can be done by using Analyse() class function. The result of 
7 // the analysis (histograms/graphs) are stored in the folder which is
8 // a data member of AliPerformanceTPC.
9 //
10 // Author: J.Otwinowski 04/02/2008 
11 // Changes by M.Knichel 15/10/2010
12 //------------------------------------------------------------------------------
13
14 /*
15  
16   // after running comparison task, read the file, and get component
17   gROOT->LoadMacro("$ALICE_ROOT/PWGPP/Macros/LoadMyLibs.C");
18   LoadMyLibs();
19
20   TFile f("Output.root");
21   AliPerformanceTPC * compObj = (AliPerformanceTPC*)coutput->FindObject("AliPerformanceTPC");
22  
23   // analyse comparison data
24   compObj->Analyse();
25
26   // the output histograms/graphs will be stored in the folder "folderTPC" 
27   compObj->GetAnalysisFolder()->ls("*");
28
29   // user can save whole comparison object (or only folder with anlysed histograms) 
30   // in the seperate output file (e.g.)
31   TFile fout("Analysed_TPC.root","recreate");
32   compObj->Write(); // compObj->GetAnalysisFolder()->Write();
33   fout.Close();
34
35 */
36
37 #include "TCanvas.h"
38 #include "TH1.h"
39 #include "TH2.h"
40 #include "TH3.h"
41 #include "TAxis.h"
42 #include "TPostScript.h"
43 #include "TString.h"
44 #include "TUUID.h"
45 #include "TTree.h"
46 #include "TChain.h"
47 #include "AliTPCPerformanceSummary.h"
48 #include "TSystem.h"
49
50 #include "AliPerformanceTPC.h" 
51 #include "AliESDEvent.h" 
52 #include "AliESDVertex.h"
53 #include "AliESDtrack.h"
54 #include "AliLog.h" 
55 #include "AliMCEvent.h" 
56 #include "AliHeader.h" 
57 #include "AliGenEventHeader.h" 
58 #include "AliStack.h" 
59 #include "AliMCInfoCuts.h" 
60 #include "AliRecInfoCuts.h" 
61 #include "AliTracker.h" 
62 #include "AliTreeDraw.h" 
63 #include "AliTPCTransform.h" 
64 #include "AliTPCseed.h" 
65 #include "AliTPCcalibDB.h" 
66 #include "AliESDfriend.h" 
67 #include "AliESDfriendTrack.h" 
68 #include "AliTPCclusterMI.h" 
69
70 using namespace std;
71
72 ClassImp(AliPerformanceTPC)
73
74 Bool_t AliPerformanceTPC::fgMergeTHnSparse = kFALSE;
75 Bool_t AliPerformanceTPC::fgUseMergeTHnSparse = kFALSE;
76
77
78 //_____________________________________________________________________________
79 /*
80 AliPerformanceTPC::AliPerformanceTPC():
81   AliPerformanceObject("AliPerformanceTPC"),
82   fTPCClustHisto(0),
83   fTPCEventHisto(0),
84   fTPCTrackHisto(0),
85   fFolderObj(0),
86
87   // Cuts 
88   fCutsRC(0),  
89   fCutsMC(0),  
90
91   // histogram folder 
92   fAnalysisFolder(0),
93   
94   fUseHLT(kFALSE)
95
96 {
97   Init();
98 }
99 */
100
101 //_____________________________________________________________________________
102 AliPerformanceTPC::AliPerformanceTPC(const Char_t* name, const Char_t* title,Int_t analysisMode,Bool_t hptGenerator, Int_t run, Bool_t highMult):
103   AliPerformanceObject(name,title,run,highMult),
104   fTPCClustHisto(0),
105   fTPCEventHisto(0),
106   fTPCTrackHisto(0),
107   fFolderObj(0),
108
109   // Cuts 
110   fCutsRC(0),  
111   fCutsMC(0),  
112
113   // histogram folder 
114   fAnalysisFolder(0),
115   
116   fUseHLT(kFALSE)
117
118 {
119   // named constructor  
120   // 
121   SetAnalysisMode(analysisMode);
122   SetHptGenerator(hptGenerator);
123
124   Init();
125 }
126
127
128 //_____________________________________________________________________________
129 AliPerformanceTPC::~AliPerformanceTPC()
130 {
131   // destructor
132    
133   if(fTPCClustHisto) delete fTPCClustHisto; fTPCClustHisto=0;     
134   if(fTPCEventHisto) delete fTPCEventHisto; fTPCEventHisto=0;     
135   if(fTPCTrackHisto) delete fTPCTrackHisto; fTPCTrackHisto=0;   
136   if(fAnalysisFolder) delete fAnalysisFolder; fAnalysisFolder=0;
137   if(fFolderObj) delete fFolderObj; fFolderObj=0;
138 }
139
140
141 //_____________________________________________________________________________
142 void AliPerformanceTPC::Init()
143 {
144   //
145   // histogram bining
146   //
147
148  
149   // set pt bins
150   Int_t nPtBins = 50;
151   Double_t ptMin = 1.e-2, ptMax = 20.;
152
153   Double_t *binsPt = 0;
154
155   if (IsHptGenerator())  { 
156         ptMax = 100.;
157   } 
158    binsPt = CreateLogAxis(nPtBins,ptMin,ptMax);
159
160
161   /*
162   const Int_t  nCOverPtBins = 80;
163   Double_t coverptMin = -10, coverptMax = 10;
164   Double_t *binsCOverPtP = 0;
165   Double_t *binsCOverPt = new Double_t[nCOverPtBins+1];
166   binsCOverPtP = CreateLogAxis(nCOverPtBins/2,0.04,coverptMax-0.04);
167   for(Int_t i=0; i < nCOverPtBins/2; i++){
168     binsCOverPt[nCOverPtBins - i] = binsCOverPtP[nCOverPtBins/2  - i];
169     binsCOverPt[i] = 0 - binsCOverPtP[nCOverPtBins/2  - i];
170  }
171  */
172
173   /*
174   Int_t nPtBins = 31;
175   Double_t binsPt[32] = {0.,0.05,0.1,0.15,0.2,0.25,0.3,0.35,0.4,0.45,0.5,0.55,0.6,0.7,0.8,0.9,1.0,1.2,1.4,1.6,1.8,2.0,2.25,2.5,2.75,3.,3.5,4.,5.,6.,8.,10.};
176   Double_t ptMin = 0., ptMax = 10.; 
177
178   if(IsHptGenerator() == kTRUE) {
179     nPtBins = 100;
180     ptMin = 0.; ptMax = 100.; 
181   }
182   */
183   //
184
185   //
186   //padRow:phi:TPCSide
187   Int_t binsTPCClustHisto[3] =   {160,  144,  2};
188   Double_t minTPCClustHisto[3] = {0.,   0.,   0.};
189   Double_t maxTPCClustHisto[3] = {160., 2.*TMath::Pi(), 2.};
190
191   fTPCClustHisto = new THnSparseF("fTPCClustHisto","padRow:phi:TPCSide",3,binsTPCClustHisto,minTPCClustHisto,maxTPCClustHisto);
192   fTPCClustHisto->GetAxis(0)->SetTitle("padRow");
193   fTPCClustHisto->GetAxis(1)->SetTitle("phi (rad)");
194   fTPCClustHisto->GetAxis(2)->SetTitle("TPCSide");
195   //fTPCClustHisto->Sumw2();
196
197   //padRow:phi:TPCSide:pad:detector:glZ
198   /*
199   Int_t binsTPCClustHisto[6] =   {160,  144,  2, 128, 72, 50};
200   Double_t minTPCClustHisto[6] = {0.,   0.,   0., 0, 0, -250};
201   Double_t maxTPCClustHisto[6] = {160., 2.*TMath::Pi(), 2., 128, 72,250};
202
203   fTPCClustHisto = new THnSparseF("fTPCClustHisto","padRow:phi:TPCSide:pad:detector:gZ",6,binsTPCClustHisto,minTPCClustHisto,maxTPCClustHisto);
204   fTPCClustHisto->GetAxis(0)->SetTitle("padRow");
205   fTPCClustHisto->GetAxis(1)->SetTitle("phi (rad)");
206   fTPCClustHisto->GetAxis(2)->SetTitle("TPCSide");
207   fTPCClustHisto->GetAxis(3)->SetTitle("pad");
208   fTPCClustHisto->GetAxis(4)->SetTitle("detector");
209   fTPCClustHisto->GetAxis(5)->SetTitle("glZ (cm)");
210   //fTPCClustHisto->Sumw2();
211   */
212   
213   Float_t scaleVxy = 1.0;
214   if(fAnalysisMode !=0) scaleVxy = 0.1; 
215
216   Int_t maxMult;
217   if (fHighMultiplicity) { maxMult = 4001; scaleVxy = 0.1;} else { maxMult = 151; }
218   // Xv:Yv:Zv:mult:multP:multN:vertStatus
219   Int_t binsTPCEventHisto[7]=  {100,  100,   100,  maxMult,  maxMult,  maxMult, 2   };
220   Double_t minTPCEventHisto[7]={-10.*scaleVxy, -10.*scaleVxy, -30., -0.5,  -0.5,  -0.5, -0.5  };
221   Double_t maxTPCEventHisto[7]={ 10.*scaleVxy,  10.*scaleVxy,  30.,  maxMult-0.5,  maxMult-0.5, maxMult-0.5, 1.5 };
222
223   fTPCEventHisto = new THnSparseF("fTPCEventHisto","Xv:Yv:Zv:mult:multP:multN:vertStatus",7,binsTPCEventHisto,minTPCEventHisto,maxTPCEventHisto);
224   fTPCEventHisto->GetAxis(0)->SetTitle("Xv (cm)");
225   fTPCEventHisto->GetAxis(1)->SetTitle("Yv (cm)");
226   fTPCEventHisto->GetAxis(2)->SetTitle("Zv (cm)");
227   fTPCEventHisto->GetAxis(3)->SetTitle("mult");
228   fTPCEventHisto->GetAxis(4)->SetTitle("multP");
229   fTPCEventHisto->GetAxis(5)->SetTitle("multN");
230   fTPCEventHisto->GetAxis(6)->SetTitle("vertStatus");
231   //fTPCEventHisto->Sumw2();
232
233   Float_t scaleDCA = 1.0;
234   if(fAnalysisMode !=0) scaleDCA = 0.1; 
235   // nTPCClust:chi2PerTPCClust:nTPCClustFindRatio:DCAr:DCAz:eta:phi:pt:charge:vertStatus
236    //Int_t binsTPCTrackHisto[10]=  { 160,  20,  60,  30, 30,  30,   144,             nPtBins,   nCOverPtBins, 2 };
237    //Double_t minTPCTrackHisto[10]={ 0.,   0.,  0., -3*scaleDCA, -3.*scaleDCA, -1.5, 0.,             ptMin,   coverptMin, -0.5 };
238    //Double_t maxTPCTrackHisto[10]={ 160., 5., 1.2, 3*scaleDCA,  3.*scaleDCA,  1.5, 2.*TMath::Pi(), ptMax,    coverptMax,  1.5 };
239    Int_t binsTPCTrackHisto[10]=  { 160,  20,  60,  30, 30,  30,   144,             nPtBins,   3, 2 };
240    Double_t minTPCTrackHisto[10]={ 0.,   0.,  0., -3*scaleDCA, -3.*scaleDCA, -1.5, 0.,             ptMin,  -1.5, -0.5 };
241    Double_t maxTPCTrackHisto[10]={ 160., 5., 1.2, 3*scaleDCA,  3.*scaleDCA,  1.5, 2.*TMath::Pi(), ptMax,    1.5,  1.5 };
242   
243   //fTPCTrackHisto = new THnSparseF("fTPCTrackHisto","nClust:chi2PerClust:nClust/nFindableClust:DCAr:DCAz:eta:phi:pt:charge/pt:vertStatus",10,binsTPCTrackHisto,minTPCTrackHisto,maxTPCTrackHisto);
244   fTPCTrackHisto = new THnSparseF("fTPCTrackHisto","nClust:chi2PerClust:nClust/nFindableClust:DCAr:DCAz:eta:phi:pt:charge:vertStatus",10,binsTPCTrackHisto,minTPCTrackHisto,maxTPCTrackHisto);
245   fTPCTrackHisto->SetBinEdges(7,binsPt);
246   //fTPCTrackHisto->SetBinEdges(8,binsCOverPt);
247   fTPCTrackHisto->GetAxis(0)->SetTitle("nClust");
248   fTPCTrackHisto->GetAxis(1)->SetTitle("chi2PerClust");
249   fTPCTrackHisto->GetAxis(2)->SetTitle("nClust/nFindableClust");
250   fTPCTrackHisto->GetAxis(3)->SetTitle("DCAr (cm)");
251   fTPCTrackHisto->GetAxis(4)->SetTitle("DCAz (cm)");
252   fTPCTrackHisto->GetAxis(5)->SetTitle("#eta");
253   fTPCTrackHisto->GetAxis(6)->SetTitle("#phi (rad)");
254   fTPCTrackHisto->GetAxis(7)->SetTitle("p_{T} (GeV/c)");
255   //fTPCTrackHisto->GetAxis(8)->SetTitle("charge/pt");
256   fTPCTrackHisto->GetAxis(8)->SetTitle("charge");
257   fTPCTrackHisto->GetAxis(9)->SetTitle("vertStatus");
258   //fTPCTrackHisto->Sumw2();
259
260   // Init cuts 
261   if(!fCutsMC) {
262     AliDebug(AliLog::kError, "ERROR: Cannot find AliMCInfoCuts object");
263   }
264   if(!fCutsRC) {
265     AliDebug(AliLog::kError, "ERROR: Cannot find AliRecInfoCuts object"); 
266   }
267
268   // init folder
269   fAnalysisFolder = CreateFolder("folderTPC","Analysis Resolution Folder");
270
271   //delete []binsCOverPt;
272   
273    // save merge status in object
274   fMergeTHnSparseObj = fgMergeTHnSparse;
275   
276 }
277
278
279 //_____________________________________________________________________________
280 void AliPerformanceTPC::ProcessTPC(AliStack* const stack, AliESDtrack *const esdTrack, AliESDEvent *const esdEvent, Bool_t vertStatus)
281 {
282 //
283 // fill TPC QA info
284 //
285   if(!esdEvent) return;
286   if(!esdTrack) return;
287
288   if(IsUseTOFBunchCrossing())
289     if(esdTrack->GetTOFBunchCrossing(esdEvent->GetMagneticField())!=0)
290       return;
291
292   if( IsUseTrackVertex() ) 
293   { 
294     // Relate TPC inner params to prim. vertex
295     const AliESDVertex *vtxESD = esdEvent->GetPrimaryVertexTracks();
296     Double_t x[3]; esdTrack->GetXYZ(x);
297     Double_t b[3]; AliTracker::GetBxByBz(x,b);
298     //    Bool_t isOK = esdTrack->RelateToVertexTPCBxByBz(vtxESD, b, kVeryBig);
299     Bool_t isOK=kFALSE;
300     if(fabs(b[2])>0.000001)
301      isOK = esdTrack->RelateToVertexTPCBxByBz(vtxESD, b, kVeryBig);
302     if(!isOK) return;
303
304     /*
305       // JMT -- recaluclate DCA for HLT if not present
306       if ( dca[0] == 0. && dca[1] == 0. ) {
307         track->GetDZ( vtxESD->GetX(), vtxESD->GetY(), vtxESD->GetZ(), esdEvent->GetMagneticField(), dca );
308       }
309     */
310   }
311
312   // Fill TPC only resolution comparison information 
313   const AliExternalTrackParam *track = esdTrack->GetTPCInnerParam();
314   if(!track) return;
315
316   Float_t dca[2], cov[3]; // dca_xy, dca_z, sigma_xy, sigma_xy_z, sigma_z
317   esdTrack->GetImpactParametersTPC(dca,cov);
318
319   Float_t q = esdTrack->Charge();
320   Float_t pt = track->Pt();
321   Float_t eta = track->Eta();
322   Float_t phi = track->Phi();
323   Int_t nClust = esdTrack->GetTPCclusters(0);
324   Int_t nFindableClust = esdTrack->GetTPCNclsF();
325
326   Float_t chi2PerCluster = 0.;
327   if(nClust>0.) chi2PerCluster = esdTrack->GetTPCchi2()/Float_t(nClust);
328
329   Float_t clustPerFindClust = 0.;
330   if(nFindableClust>0.) clustPerFindClust = Float_t(nClust)/nFindableClust;
331   
332   Float_t qpt = 0;
333   if( fabs(pt)>0 ) qpt = q/fabs(pt);
334
335   // filter out noise tracks
336   if(esdTrack->GetTPCsignal() < 5) return;
337
338   //
339   // select primaries
340   //
341   Double_t dcaToVertex = -1;
342   if( fCutsRC->GetDCAToVertex2D() ) 
343   {
344       dcaToVertex = TMath::Sqrt(dca[0]*dca[0]/fCutsRC->GetMaxDCAToVertexXY()/fCutsRC->GetMaxDCAToVertexXY() + dca[1]*dca[1]/fCutsRC->GetMaxDCAToVertexZ()/fCutsRC->GetMaxDCAToVertexZ()); 
345   }
346   if(fCutsRC->GetDCAToVertex2D() && dcaToVertex > 1) return;
347   if(!fCutsRC->GetDCAToVertex2D() && TMath::Abs(dca[0]) > fCutsRC->GetMaxDCAToVertexXY()) return;
348   if(!fCutsRC->GetDCAToVertex2D() && TMath::Abs(dca[1]) > fCutsRC->GetMaxDCAToVertexZ()) return;
349
350   //Double_t vTPCTrackHisto[10] = {nClust,chi2PerCluster,clustPerFindClust,dca[0],dca[1],eta,phi,pt,qpt,vertStatus};
351   Double_t vTPCTrackHisto[10] = {static_cast<Double_t>(nClust),static_cast<Double_t>(chi2PerCluster),static_cast<Double_t>(clustPerFindClust),static_cast<Double_t>(dca[0]),static_cast<Double_t>(dca[1]),static_cast<Double_t>(eta),static_cast<Double_t>(phi),static_cast<Double_t>(pt),static_cast<Double_t>(q),static_cast<Double_t>(vertStatus)};
352   fTPCTrackHisto->Fill(vTPCTrackHisto); 
353  
354   //
355   // Fill rec vs MC information
356   //
357   if(!stack) return;
358
359 }
360
361
362 //_____________________________________________________________________________
363 void AliPerformanceTPC::ProcessTPCITS(AliStack* const stack, AliESDtrack *const esdTrack, AliESDEvent* const esdEvent, Bool_t vertStatus)
364 {
365   // Fill comparison information (TPC+ITS) 
366   if(!esdTrack) return;
367   if(!esdEvent) return;
368
369   if( IsUseTrackVertex() ) 
370   { 
371     // Relate TPC inner params to prim. vertex
372     const AliESDVertex *vtxESD = esdEvent->GetPrimaryVertexTracks();
373     Double_t x[3]; esdTrack->GetXYZ(x);
374     Double_t b[3]; AliTracker::GetBxByBz(x,b);
375     Bool_t isOK = esdTrack->RelateToVertexBxByBz(vtxESD, b, kVeryBig);
376     if(!isOK) return;
377
378     /*
379       // JMT -- recaluclate DCA for HLT if not present
380       if ( dca[0] == 0. && dca[1] == 0. ) {
381         track->GetDZ( vtxESD->GetX(), vtxESD->GetY(), vtxESD->GetZ(), esdEvent->GetMagneticField(), dca );
382       }
383     */
384   }
385
386   Float_t dca[2], cov[3]; // dca_xy, dca_z, sigma_xy, sigma_xy_z, sigma_z
387   esdTrack->GetImpactParameters(dca,cov);
388
389   if ((esdTrack->GetStatus()&AliESDtrack::kITSrefit)==0) return; // ITS refit
390   if ((esdTrack->GetStatus()&AliESDtrack::kTPCrefit)==0) return; // TPC refit
391   if ((esdTrack->HasPointOnITSLayer(0)==kFALSE)&&(esdTrack->HasPointOnITSLayer(1)==kFALSE)) return; // at least one SPD
392   //if (esdTrack->GetITSclusters(0)<fCutsRC->GetMinNClustersITS()) return;  // min. nb. ITS clusters
393
394   Float_t q = esdTrack->Charge();
395   Float_t pt = esdTrack->Pt();
396   Float_t eta = esdTrack->Eta();
397   Float_t phi = esdTrack->Phi();
398   Int_t nClust = esdTrack->GetTPCclusters(0);
399   Int_t nFindableClust = esdTrack->GetTPCNclsF();
400
401   Float_t chi2PerCluster = 0.;
402   if(nClust>0.) chi2PerCluster = esdTrack->GetTPCchi2()/Float_t(nClust);
403
404   Float_t clustPerFindClust = 0.;
405   if(nFindableClust>0.) clustPerFindClust = Float_t(nClust)/nFindableClust;
406   
407   Float_t qpt = 0;
408   if( fabs(pt)>0 ) qpt = q/fabs(pt);
409
410   //
411   // select primaries
412   //
413   Double_t dcaToVertex = -1;
414   if( fCutsRC->GetDCAToVertex2D() ) 
415   {
416       dcaToVertex = TMath::Sqrt(dca[0]*dca[0]/fCutsRC->GetMaxDCAToVertexXY()/fCutsRC->GetMaxDCAToVertexXY() + dca[1]*dca[1]/fCutsRC->GetMaxDCAToVertexZ()/fCutsRC->GetMaxDCAToVertexZ()); 
417   }
418   if(fCutsRC->GetDCAToVertex2D() && dcaToVertex > 1) return;
419   if(!fCutsRC->GetDCAToVertex2D() && TMath::Abs(dca[0]) > fCutsRC->GetMaxDCAToVertexXY()) return;
420   if(!fCutsRC->GetDCAToVertex2D() && TMath::Abs(dca[1]) > fCutsRC->GetMaxDCAToVertexZ()) return;
421
422   Double_t vTPCTrackHisto[10] = {static_cast<Double_t>(nClust),static_cast<Double_t>(chi2PerCluster),static_cast<Double_t>(clustPerFindClust),static_cast<Double_t>(dca[0]),static_cast<Double_t>(dca[1]),static_cast<Double_t>(eta),static_cast<Double_t>(phi),static_cast<Double_t>(pt),static_cast<Double_t>(q),static_cast<Double_t>(vertStatus)};
423   fTPCTrackHisto->Fill(vTPCTrackHisto); 
424  
425   //
426   // Fill rec vs MC information
427   //
428   if(!stack) return;
429 }
430
431
432 //_____________________________________________________________________________
433 void AliPerformanceTPC::ProcessConstrained(AliStack* const /*stack*/, AliESDtrack *const /*esdTrack*/, AliESDEvent* const /*esdEvent*/)
434 {
435   // Fill comparison information (constarained parameters) 
436   AliDebug(AliLog::kWarning, "Warning: Not implemented");
437 }
438
439
440 //_____________________________________________________________________________
441 void AliPerformanceTPC::Exec(AliMCEvent* const mcEvent, AliESDEvent *const esdEvent, AliESDfriend *const esdFriend, const Bool_t bUseMC, const Bool_t bUseESDfriend)
442 {
443   // Process comparison information 
444   //
445
446   if(!esdEvent) 
447   {
448     Error("Exec","esdEvent not available");
449     return;
450   }
451   AliHeader* header = 0;
452   AliGenEventHeader* genHeader = 0;
453   AliStack* stack = 0;
454   TArrayF vtxMC(3);
455   
456   if(bUseMC)
457   {
458     if(!mcEvent) {
459       Error("Exec","mcEvent not available");
460       return;
461     }
462     // get MC event header
463     header = mcEvent->Header();
464     if (!header) {
465       Error("Exec","Header not available");
466       return;
467     }
468     // MC particle stack
469     stack = mcEvent->Stack();
470     if (!stack) {
471       Error("Exec","Stack not available");
472       return;
473     }
474     // get MC vertex
475     genHeader = header->GenEventHeader();
476     if (!genHeader) {
477       Error("Exec","Could not retrieve genHeader from Header");
478       return;
479     }
480     genHeader->PrimaryVertex(vtxMC);
481   } 
482   
483   // trigger
484   if(!bUseMC &&GetTriggerClass()) {
485     Bool_t isEventTriggered = esdEvent->IsTriggerClassFired(GetTriggerClass());
486     if(!isEventTriggered) return; 
487   }
488
489   // get TPC event vertex
490   const AliESDVertex *vtxESD = NULL; 
491   if(fUseTrackVertex) {
492     vtxESD = esdEvent->GetPrimaryVertexTracks();
493   } else {
494     vtxESD = esdEvent->GetPrimaryVertexTPC();
495   }
496   if(!vtxESD) return;
497
498   //  events with rec. vertex
499   Int_t mult=0; Int_t multP=0; Int_t multN=0;
500   
501   // changed to take all events but store vertex status
502 //  if(vtxESD->GetStatus() >0)
503 //  {
504   // store vertex status
505   Bool_t vertStatus = vtxESD->GetStatus();
506   //  Process ESD events
507   for (Int_t iTrack = 0; iTrack < esdEvent->GetNumberOfTracks(); iTrack++) 
508   { 
509     AliESDtrack *track = esdEvent->GetTrack(iTrack);
510     if(!track) continue;
511
512     // if not fUseKinkDaughters don't use tracks with kink index > 0
513     if(!fUseKinkDaughters && track->GetKinkIndex(0) > 0) continue;
514     
515     if(bUseESDfriend && esdFriend && esdFriend->TestSkipBit()==kFALSE && iTrack<esdFriend->GetNumberOfTracks()) 
516     {
517       AliESDfriendTrack *friendTrack=esdFriend->GetTrack(iTrack);
518       if(friendTrack) 
519       {
520         //
521         TObject *calibObject=0;
522         AliTPCseed *seed=0;
523         for (Int_t j=0;(calibObject=friendTrack->GetCalibObject(j));++j) {
524             if ((seed=dynamic_cast<AliTPCseed*>(calibObject))) {
525             break;
526           }
527         }
528
529         // 
530         for (Int_t irow=0;irow<159;irow++) {
531         if(!seed) continue;
532           
533           AliTPCclusterMI *cluster=seed->GetClusterPointer(irow);
534           if (!cluster) continue;
535
536              Float_t gclf[3];
537              cluster->GetGlobalXYZ(gclf);
538
539              //Double_t x[3]={cluster->GetRow(),cluster->GetPad(),cluster->GetTimeBin()};
540              //Int_t i[1]={cluster->GetDetector()};
541              //transform->Transform(x,i,0,1);
542              //printf("gx %f gy  %f  gz %f \n", cluster->GetX(), cluster->GetY(),cluster->GetZ());
543              //printf("gclf[0] %f gclf[1]  %f  gclf[2] %f \n", gclf[0], gclf[1],  gclf[2]);
544      
545              Int_t TPCside; 
546              if(gclf[2]>0.) TPCside=0; // A side 
547              else TPCside=1;
548
549              //
550              //Double_t vTPCClust1[3] = { gclf[0], gclf[1],  TPCside };
551              //fTPCClustHisto1->Fill(vTPCClust1);
552
553              //  
554              Double_t phi = TMath::ATan2(gclf[1],gclf[0]);
555              if(phi < 0) phi += 2.*TMath::Pi();
556             
557              //Float_t pad = cluster->GetPad();
558              //Int_t detector = cluster->GetDetector();
559              //Double_t vTPCClust[6] = { irow, phi, TPCside, pad, detector, gclf[2] };
560              Double_t vTPCClust[3] = { static_cast<Double_t>(irow), phi, static_cast<Double_t>(TPCside) };
561              fTPCClustHisto->Fill(vTPCClust);
562         }
563       }
564     }
565
566     if(GetAnalysisMode() == 0) ProcessTPC(stack,track,esdEvent,vertStatus);
567     else if(GetAnalysisMode() == 1) ProcessTPCITS(stack,track,esdEvent,vertStatus);
568     else if(GetAnalysisMode() == 2) ProcessConstrained(stack,track,esdEvent);
569     else {
570       printf("ERROR: AnalysisMode %d \n",fAnalysisMode);
571       return;
572     }
573
574     // TPC only
575     if(!fUseHLT){
576       if(GetAnalysisMode() == 0) {
577         AliESDtrack *tpcTrack = AliESDtrackCuts::GetTPCOnlyTrack(esdEvent,iTrack);
578         if(!tpcTrack) continue;
579         // track selection
580         if( fCutsRC->AcceptTrack(tpcTrack) ) { 
581           mult++;
582           if(tpcTrack->Charge()>0.) multP++;
583           if(tpcTrack->Charge()<0.) multN++;
584         }
585         if(tpcTrack) delete tpcTrack;
586       }
587       else {
588        // track selection
589         if( fCutsRC->AcceptTrack(track) ) { 
590           mult++;
591           if(track->Charge()>0.) multP++;
592           if(track->Charge()<0.) multN++;
593         }
594       } 
595     }
596     else {
597       if( fCutsRC->AcceptTrack(track) ) { 
598         //Printf("Still here for HLT");
599         mult++;
600         if(track->Charge()>0.) multP++;
601         if(track->Charge()<0.) multN++;
602       }
603     }
604   }
605
606   Double_t vTPCEvent[7] = {vtxESD->GetXv(),vtxESD->GetYv(),vtxESD->GetZv(),static_cast<Double_t>(mult),static_cast<Double_t>(multP),static_cast<Double_t>(multN),static_cast<Double_t>(vtxESD->GetStatus())};
607   fTPCEventHisto->Fill(vTPCEvent);
608 }
609
610
611 //_____________________________________________________________________________
612 void AliPerformanceTPC::Analyse()
613 {
614     //
615     // Analyse comparison information and store output histograms
616     // in the folder "folderTPC"
617     //
618     TH1::AddDirectory(kFALSE);
619     TH1::SetDefaultSumw2(kFALSE);
620     TObjArray *aFolderObj = new TObjArray;
621     //aFolderObj->SetOwner(); // objects are owned by fanalysisFolder
622     TString selString;
623
624     //
625     // Cluster histograms
626     //
627     AddProjection(aFolderObj, "clust", fTPCClustHisto, 0, 1, 2);
628     
629     selString = "all";
630     for(Int_t i=0; i <= 2; i++) {
631         AddProjection(aFolderObj, "clust", fTPCClustHisto, i, &selString);
632     }
633     
634     //fTPCClustHisto->GetAxis(2)->SetRange(1,1); // A-side
635     //selString = "A_side";
636     //AddProjection(aFolderObj, fTPCClustHisto, 0, 1, &selString);
637     
638     //fTPCClustHisto->GetAxis(2)->SetRange(2,2); // C-side
639     //selString = "C_side";
640     //AddProjection(aFolderObj, fTPCClustHisto, 0, 1, &selString);
641     
642     //reset range
643     fTPCClustHisto->GetAxis(2)->SetRange(1,2); 
644     
645     //
646     // event histograms
647     //
648     for(Int_t i=0; i<=6; i++) {
649       AddProjection(aFolderObj, "event", fTPCEventHisto, i);
650     }    
651     AddProjection(aFolderObj, "event", fTPCEventHisto, 4, 5);
652     AddProjection(aFolderObj, "event", fTPCEventHisto, 0, 1);
653     AddProjection(aFolderObj, "event", fTPCEventHisto, 0, 3);
654     AddProjection(aFolderObj, "event", fTPCEventHisto, 1, 3);
655     AddProjection(aFolderObj, "event", fTPCEventHisto, 2, 3);
656
657     // reconstructed vertex status > 0
658     fTPCEventHisto->GetAxis(6)->SetRange(2,2);
659     selString = "recVertex";
660     for(Int_t i=0; i<=5; i++) {
661       AddProjection(aFolderObj, "event", fTPCEventHisto, i, &selString);
662     }
663     AddProjection(aFolderObj, "event", fTPCEventHisto, 4, 5, &selString);
664     AddProjection(aFolderObj, "event", fTPCEventHisto, 0, 1, &selString);
665     AddProjection(aFolderObj, "event", fTPCEventHisto, 0, 3, &selString);
666     AddProjection(aFolderObj, "event", fTPCEventHisto, 1, 3, &selString);
667     AddProjection(aFolderObj, "event", fTPCEventHisto, 2, 3, &selString);
668
669     // reset cuts
670     fTPCEventHisto->GetAxis(6)->SetRange(1,2);
671
672     //
673     // Track histograms 
674     // 
675     // all with vertex
676     fTPCTrackHisto->GetAxis(8)->SetRangeUser(-1.5,1.5);
677     fTPCTrackHisto->GetAxis(9)->SetRangeUser(0.5,1.5);
678     selString = "all_recVertex";
679     for(Int_t i=0; i <= 9; i++) {
680         AddProjection(aFolderObj, "track", fTPCTrackHisto, i, &selString);        
681     }
682
683      AddProjection(aFolderObj, "track", fTPCTrackHisto, 5, 8, &selString); 
684
685     for(Int_t i=0; i <= 4; i++) {
686         AddProjection(aFolderObj, "track", fTPCTrackHisto, i, 5, 7, &selString);        
687     }    
688
689
690
691     // Track histograms (pos with vertex)
692     fTPCTrackHisto->GetAxis(8)->SetRangeUser(0,1.5);
693     selString = "pos_recVertex";
694     for(Int_t i=0; i <= 9; i++) {
695         AddProjection(aFolderObj, "track", fTPCTrackHisto, i, &selString);
696     }
697     for(Int_t i=0; i <= 4; i++) { for(Int_t j=5; j <= 5; j++) { for(Int_t k=j+1; k <= 7; k++) {
698         AddProjection(aFolderObj, "track", fTPCTrackHisto, i, j, k, &selString);
699     }  }  }
700     AddProjection(aFolderObj, "track", fTPCTrackHisto, 0, 1, 2, &selString);
701     AddProjection(aFolderObj, "track", fTPCTrackHisto, 0, 1, 5, &selString);
702     AddProjection(aFolderObj, "track", fTPCTrackHisto, 0, 2, 5, &selString);
703     AddProjection(aFolderObj, "track", fTPCTrackHisto, 1, 2, 5, &selString);
704     AddProjection(aFolderObj, "track", fTPCTrackHisto, 3, 4, 5, &selString);
705     AddProjection(aFolderObj, "track", fTPCTrackHisto, 5, 6, 7, &selString);
706   
707     // Track histograms (neg with vertex)
708     fTPCTrackHisto->GetAxis(8)->SetRangeUser(-1.5,0);
709     selString = "neg_recVertex";
710     for(Int_t i=0; i <= 9; i++) {
711         AddProjection(aFolderObj, "track", fTPCTrackHisto, i, &selString);
712     }
713     for(Int_t i=0; i <= 4; i++) { for(Int_t j=5; j <= 5; j++) { for(Int_t k=j+1; k <= 7; k++) {
714         AddProjection(aFolderObj, "track", fTPCTrackHisto, i, j, k, &selString);
715     }  }  }
716     AddProjection(aFolderObj, "track", fTPCTrackHisto, 0, 1, 2, &selString);
717     AddProjection(aFolderObj, "track", fTPCTrackHisto, 0, 1, 5, &selString);
718     AddProjection(aFolderObj, "track", fTPCTrackHisto, 0, 2, 5, &selString);
719     AddProjection(aFolderObj, "track", fTPCTrackHisto, 1, 2, 5, &selString);
720     AddProjection(aFolderObj, "track", fTPCTrackHisto, 3, 4, 5, &selString);
721     AddProjection(aFolderObj, "track", fTPCTrackHisto, 5, 6, 7, &selString);
722
723     //restore cuts
724     fTPCTrackHisto->GetAxis(8)->SetRangeUser(-1.5,1.5);
725     fTPCTrackHisto->GetAxis(9)->SetRangeUser(-0.5,1.5);
726   
727   
728     printf("exportToFolder\n");
729     // export objects to analysis folder
730     fAnalysisFolder = ExportToFolder(aFolderObj);
731     if (fFolderObj) delete fFolderObj;
732     fFolderObj = aFolderObj;
733     aFolderObj=0;
734 }
735
736
737 //_____________________________________________________________________________
738 TFolder* AliPerformanceTPC::ExportToFolder(TObjArray * array) 
739 {
740   // recreate folder avery time and export objects to new one
741   //
742   AliPerformanceTPC * comp=this;
743   TFolder *folder = comp->GetAnalysisFolder();
744
745   TString name, title;
746   TFolder *newFolder = 0;
747   Int_t i = 0;
748   Int_t size = array->GetSize();
749
750   if(folder) { 
751      // get name and title from old folder
752      name = folder->GetName();  
753      title = folder->GetTitle();  
754
755          // delete old one
756      delete folder;
757
758          // create new one
759      newFolder = CreateFolder(name.Data(),title.Data());
760      newFolder->SetOwner();
761
762          // add objects to folder
763      while(i < size) {
764            newFolder->Add(array->At(i));
765            i++;
766          }
767   }
768
769 return newFolder;
770 }
771
772 //_____________________________________________________________________________
773 Long64_t AliPerformanceTPC::Merge(TCollection* const list) 
774 {
775   // Merge list of objects (needed by PROOF)
776
777   if (!list)
778   return 0;
779
780   if (list->IsEmpty())
781   return 1;
782   
783   Bool_t merge = ((fgUseMergeTHnSparse && fgMergeTHnSparse) || (!fgUseMergeTHnSparse && fMergeTHnSparseObj));
784
785   TIterator* iter = list->MakeIterator();
786   TObject* obj = 0;
787   TObjArray* objArrayList = 0;
788   objArrayList = new TObjArray();
789
790   // collection of generated histograms
791   Int_t count=0;
792   while((obj = iter->Next()) != 0) 
793   {
794     AliPerformanceTPC* entry = dynamic_cast<AliPerformanceTPC*>(obj);
795     if (entry == 0) continue; 
796     if (merge) {
797         if ((fTPCClustHisto) && (entry->fTPCClustHisto)) { fTPCClustHisto->Add(entry->fTPCClustHisto); }
798         if ((fTPCEventHisto) && (entry->fTPCEventHisto)) { fTPCEventHisto->Add(entry->fTPCEventHisto); }
799         if ((fTPCTrackHisto) && (entry->fTPCTrackHisto)) { fTPCTrackHisto->Add(entry->fTPCTrackHisto); }
800     }
801     // the analysisfolder is only merged if present
802     if (entry->fFolderObj) { objArrayList->Add(entry->fFolderObj); }
803
804     count++;
805   }
806   if (fFolderObj) { fFolderObj->Merge(objArrayList); } 
807   // to signal that track histos were not merged: reset
808   if (!merge) { fTPCTrackHisto->Reset(); fTPCClustHisto->Reset(); fTPCEventHisto->Reset(); }
809   // delete
810   if (objArrayList)  delete objArrayList;  objArrayList=0;
811 return count;
812 }
813
814
815 //_____________________________________________________________________________
816 TFolder* AliPerformanceTPC::CreateFolder(TString name, TString title) 
817
818 // create folder for analysed histograms
819 //
820 TFolder *folder = 0;
821   folder = new TFolder(name.Data(),title.Data());
822
823   return folder;
824 }
825
826 //_____________________________________________________________________________
827 TTree* AliPerformanceTPC::CreateSummary()
828 {
829     // implementaion removed, switched back to use AliPerformanceSummary (now called in AliPerformanceTask)
830     return 0;
831 }
832