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