fd6442cbac668edf019f1f9f05f9e49333dd59ee
[u/mrichter/AliRoot.git] / TRD / qaRec / AliTRDcheckDET.cxx
1 ////////////////////////////////////////////////////////////////////////////
2 //                                                                        //
3 //                                                                        //
4 //  Basic checks for tracking and detector performance                    //
5 //  
6 //     For the moment (15.10.2009) the following checks are implemented    //
7 //       - Number of clusters/track
8 //       - Number of clusters/tracklet
9 //       - Number of tracklets/track from different sources (Barrel, StandAlone)
10 //       - Number of findable tracklets
11 //       - Number of tracks per event and TRD sector
12 //       - <PH>
13 //       - Chi2 distribution for tracks
14 //       - Charge distribution per cluster and tracklet
15 //       - Number of events and tracks per trigger source 
16 //       - Trigger purity
17 //       - Track and Tracklet propagation status
18 //
19 //  Authors:                                                              //
20 //    Anton Andronic <A.Andronic@gsi.de>                                  //
21 //    Alexandru Bercuci <A.Bercuci@gsi.de>                                //
22 //    Markus Fasel <M.Fasel@gsi.de>                                       //
23 //                                                                        //
24 ////////////////////////////////////////////////////////////////////////////
25
26 #include <TAxis.h>
27 #include <TCanvas.h>
28 #include <TFile.h>
29 #include <TH1F.h>
30 #include <TH1I.h>
31 #include <TF1.h>
32 #include <TGaxis.h>
33 #include <TGraph.h>
34 #include <TLegend.h>
35 #include <TMath.h>
36 #include <TMap.h>
37 #include <TObjArray.h>
38 #include <TObject.h>
39 #include <TObjString.h>
40
41 #include <TPad.h>
42 #include <TProfile.h>
43 #include <TProfile2D.h>
44 #include <TROOT.h>
45
46 #include "AliLog.h"
47 #include "AliTRDcluster.h"
48 #include "AliESDHeader.h"
49 #include "AliESDRun.h"
50 #include "AliESDtrack.h"
51 #include "AliTRDgeometry.h"
52 #include "AliTRDpadPlane.h"
53 #include "AliTRDSimParam.h"
54 #include "AliTRDseedV1.h"
55 #include "AliTRDtrackV1.h"
56 #include "AliTRDtrackerV1.h"
57 #include "AliTRDReconstructor.h"
58 #include "AliTrackReference.h"
59 #include "AliTrackPointArray.h"
60 #include "AliTracker.h"
61 #include "TTreeStream.h"
62
63 #include "info/AliTRDtrackInfo.h"
64 #include "info/AliTRDeventInfo.h"
65 #include "AliTRDcheckDET.h"
66
67 #include <cstdio>
68 #include <iostream>
69
70 ClassImp(AliTRDcheckDET)
71
72 //_______________________________________________________
73 AliTRDcheckDET::AliTRDcheckDET():
74   AliTRDrecoTask("checkDET", "Basic TRD data checker")
75   ,fEventInfo(0x0)
76   ,fTriggerNames(0x0)
77   ,fReconstructor(0x0)
78   ,fGeo(0x0)
79 {
80   //
81   // Default constructor
82   //
83   DefineInput(1,AliTRDeventInfo::Class());
84   fReconstructor = new AliTRDReconstructor;
85   fReconstructor->SetRecoParam(AliTRDrecoParam::GetLowFluxParam());
86   fGeo = new AliTRDgeometry;
87   InitFunctorList();
88 }
89
90 //_______________________________________________________
91 AliTRDcheckDET::~AliTRDcheckDET(){
92   //
93   // Destructor
94   // 
95   if(fTriggerNames) delete fTriggerNames;
96   delete fReconstructor;
97   delete fGeo;
98 }
99
100 //_______________________________________________________
101 void AliTRDcheckDET::ConnectInputData(Option_t *opt){
102   //
103   // Connect the Input data with the task
104   //
105   AliTRDrecoTask::ConnectInputData(opt);
106   fEventInfo = dynamic_cast<AliTRDeventInfo *>(GetInputData(1));
107 }
108
109 //_______________________________________________________
110 void AliTRDcheckDET::CreateOutputObjects(){
111   //
112   // Create Output Objects
113   //
114   OpenFile(0,"RECREATE");
115   fContainer = Histos();
116   if(!fTriggerNames) fTriggerNames = new TMap();
117 }
118
119 //_______________________________________________________
120 void AliTRDcheckDET::Exec(Option_t *opt){
121   //
122   // Execution function
123   // Filling TRD quality histos
124   //
125   if(!HasMCdata() && fEventInfo->GetEventHeader()->GetEventType() != 7) return; // For real data we select only physical events
126   AliTRDrecoTask::Exec(opt);  
127   Int_t nTracks = 0;            // Count the number of tracks per event
128   Int_t triggermask = fEventInfo->GetEventHeader()->GetTriggerMask();
129   TString triggername =  fEventInfo->GetRunInfo()->GetFiredTriggerClasses(triggermask);
130   if(fDebugLevel > 6)printf("Trigger cluster: %d, Trigger class: %s\n", triggermask, triggername.Data());
131   dynamic_cast<TH1F *>(fContainer->UncheckedAt(kNeventsTrigger))->Fill(triggermask);
132   for(Int_t iti = 0; iti < fTracks->GetEntriesFast(); iti++){
133     if(!fTracks->UncheckedAt(iti)) continue;
134     AliTRDtrackInfo *fTrackInfo = dynamic_cast<AliTRDtrackInfo *>(fTracks->UncheckedAt(iti));
135     if(!fTrackInfo->GetTrack()) continue;
136     nTracks++;
137   }
138   if(nTracks){
139     dynamic_cast<TH1F *>(fContainer->UncheckedAt(kNeventsTriggerTracks))->Fill(triggermask);
140     dynamic_cast<TH1F *>(fContainer->UncheckedAt(kNtracksEvent))->Fill(nTracks);
141   }
142   if(triggermask <= 20 && !fTriggerNames->FindObject(Form("%d", triggermask))){
143     fTriggerNames->Add(new TObjString(Form("%d", triggermask)), new TObjString(triggername));
144     // also set the label for both histograms
145     TH1 *histo = dynamic_cast<TH1F *>(fContainer->UncheckedAt(kNeventsTriggerTracks));
146     histo->GetXaxis()->SetBinLabel(histo->FindBin(triggermask), triggername);
147     histo = dynamic_cast<TH1F *>(fContainer->UncheckedAt(kNeventsTrigger));
148     histo->GetXaxis()->SetBinLabel(histo->FindBin(triggermask), triggername);
149   }
150   PostData(0, fContainer);
151 }
152
153
154 //_______________________________________________________
155 Bool_t AliTRDcheckDET::PostProcess(){
156   //
157   // Do Postprocessing (for the moment set the number of Reference histograms)
158   //
159   
160   TH1 * h = 0x0;
161   
162   // Calculate of the trigger clusters purity
163   h = dynamic_cast<TH1F *>(fContainer->FindObject("hEventsTrigger"));
164   TH1F *h1 = dynamic_cast<TH1F *>(fContainer->FindObject("hEventsTriggerTracks"));
165   h1->Divide(h);
166   Float_t purities[20], val = 0;
167   TString triggernames[20];
168   Int_t nTriggerClasses = 0;
169   for(Int_t ibin = 1; ibin <= h->GetNbinsX(); ibin++){
170     if((val = h1->GetBinContent(ibin))){
171       purities[nTriggerClasses] = val;
172       triggernames[nTriggerClasses] = h1->GetXaxis()->GetBinLabel(ibin);
173       nTriggerClasses++;
174     }
175   }
176   h = dynamic_cast<TH1F *>(fContainer->UncheckedAt(kTriggerPurity));
177   TAxis *ax = h->GetXaxis();
178   for(Int_t itrg = 0; itrg < nTriggerClasses; itrg++){
179     h->Fill(itrg, purities[itrg]);
180     ax->SetBinLabel(itrg+1, triggernames[itrg].Data());
181   }
182   ax->SetRangeUser(-0.5, nTriggerClasses+.5);
183   h->GetYaxis()->SetRangeUser(0,1);
184
185   // track status
186   h=dynamic_cast<TH1F*>(fContainer->At(kTrackStatus));
187   Float_t ok = h->GetBinContent(1);
188   Int_t nerr = h->GetNbinsX();
189   for(Int_t ierr=nerr; ierr--;){
190     h->SetBinContent(ierr+1, 1.e2*h->GetBinContent(ierr+1)/ok);
191   }
192   h->SetBinContent(1, 0.);
193
194   // tracklet status
195   TObjArray *arr = dynamic_cast<TObjArray*>(fContainer->UncheckedAt(kTrackletStatus));
196   for(Int_t ily = AliTRDgeometry::kNlayer; ily--;){
197     h=dynamic_cast<TH1F*>(arr->At(ily));
198     Float_t okB = h->GetBinContent(1);
199     Int_t nerrB = h->GetNbinsX();
200     for(Int_t ierr=nerrB; ierr--;){
201       h->SetBinContent(ierr+1, 1.e2*h->GetBinContent(ierr+1)/okB);
202     }
203     h->SetBinContent(1, 0.);
204   }
205
206   fNRefFigures = 18;
207
208   return kTRUE;
209 }
210
211 //_______________________________________________________
212 Bool_t AliTRDcheckDET::GetRefFigure(Int_t ifig){
213   //
214   // Setting Reference Figures
215   //
216   gPad->SetLogy(0);
217   gPad->SetLogx(0);
218   TH1 *h = 0x0; TObjArray *arr=0x0;
219   TLegend *leg = 0x0;
220   Bool_t kFIRST(1);
221   switch(ifig){
222   case kNclustersTrack:
223     (h=(TH1F*)fContainer->FindObject("hNcls"))->Draw("pl");
224     PutTrendValue("NClustersTrack", h->GetMean());
225     PutTrendValue("NClustersTrackRMS", h->GetRMS());
226     return kTRUE;
227   case kNclustersTracklet:
228     (h =(TH1F*)fContainer->FindObject("hNclTls"))->Draw("pc");
229     PutTrendValue("NClustersTracklet", h->GetMean());
230     PutTrendValue("NClustersTrackletRMS", h->GetRMS());
231     return kTRUE;
232   case kNtrackletsTrack:
233     h=MakePlotNTracklets();
234     PutTrendValue("NTrackletsTrack", h->GetMean());
235     PutTrendValue("NTrackletsTrackRMS", h->GetRMS());
236     return kTRUE;
237   case kNtrackletsCross:
238     h = (TH1F*)fContainer->FindObject("hNtlsCross");
239     if(!MakeBarPlot(h, kRed)) break;
240     PutTrendValue("NTrackletsCross", h->GetMean());
241     PutTrendValue("NTrackletsCrossRMS", h->GetRMS());
242     return kTRUE;
243   case kNtrackletsFindable:
244     h = (TH1F*)fContainer->FindObject("hNtlsFindable");
245     if(!MakeBarPlot(h, kGreen)) break;
246     PutTrendValue("NTrackletsFindable", h->GetMean());
247     PutTrendValue("NTrackletsFindableRMS", h->GetRMS());
248     return kTRUE;
249   case kNtracksEvent:
250     (h = (TH1F*)fContainer->FindObject("hNtrks"))->Draw("pl");
251     PutTrendValue("NTracksEvent", h->GetMean());
252     PutTrendValue("NTracksEventRMS", h->GetRMS());
253     return kTRUE;
254   case kNtracksSector:
255     h = (TH1F*)fContainer->FindObject("hNtrksSector");
256     if(!MakeBarPlot(h, kGreen)) break;
257     PutTrendValue("NTracksSector", h->Integral()/h->GetNbinsX());
258     return kTRUE;
259   case kTrackStatus:
260     if(!(h=(TH1F *)fContainer->FindObject("hTrackStatus"))) break;
261     h->GetXaxis()->SetRangeUser(0.5, -1);
262     h->GetYaxis()->CenterTitle();
263     h->Draw("c");
264     PutTrendValue("TrackStatus", h->Integral());
265     gPad->SetLogy(0);
266     return kTRUE;
267   case kTrackletStatus:
268     if(!(arr = dynamic_cast<TObjArray*>(fContainer->At(kTrackletStatus)))) break;
269     leg = new TLegend(.68, .7, .98, .98);
270     leg->SetBorderSize(1);leg->SetFillColor(0);
271     leg->SetHeader("TRD layer");
272     for(Int_t ily=0; ily<AliTRDgeometry::kNlayer; ily++){
273       if(!(h=dynamic_cast<TH1F*>(arr->At(ily)))) continue;
274       if(kFIRST){
275         h->Draw("c");
276         h->GetXaxis()->SetRangeUser(0.5, -1);
277         h->GetYaxis()->CenterTitle();
278         kFIRST = kFALSE;
279       } else h->Draw("samec");
280       leg->AddEntry(h, Form("%d", ily), "l");
281       PutTrendValue(Form("TrackletStatus%d", ily), h->Integral());
282     }
283     leg->Draw();
284     gPad->SetLogy(0);
285     return kTRUE;
286   case kChi2:
287     MakePlotChi2();
288     return kTRUE;
289   case kPH:
290     MakePlotPulseHeight();
291     gPad->SetLogy(0);
292     return kTRUE;
293   case kChargeCluster:
294     (h = (TH1F*)fContainer->FindObject("hQcl"))->Draw("c");
295     gPad->SetLogy(1);
296     PutTrendValue("ChargeCluster", h->GetMaximumBin());
297     PutTrendValue("ChargeClusterRMS", h->GetRMS());
298     return kTRUE;
299   case kChargeTracklet:
300     (h=(TH1F*)fContainer->FindObject("hQtrklt"))->Draw("c");
301     PutTrendValue("ChargeTracklet", h->GetMaximumBin());
302     PutTrendValue("ChargeTrackletRMS", h->GetRMS());
303     return kTRUE;
304   case kNeventsTrigger:
305     ((TH1F*)fContainer->FindObject("hEventsTrigger"))->Draw("");
306     return kTRUE;
307   case kNeventsTriggerTracks:
308     ((TH1F*)fContainer->FindObject("hEventsTriggerTracks"))->Draw("");
309     return kTRUE;
310   case kTriggerPurity: 
311     if(!MakeBarPlot((TH1F*)fContainer->FindObject("hTriggerPurity"), kGreen)) break;
312     break;
313   default:
314     break;
315   }
316   AliInfo(Form("Reference plot [%d] missing result", ifig));
317   return kFALSE;
318 }
319
320 //_______________________________________________________
321 TObjArray *AliTRDcheckDET::Histos(){
322   //
323   // Create QA histograms
324   //
325   if(fContainer) return fContainer;
326   
327   fContainer = new TObjArray(20);
328   //fContainer->SetOwner(kTRUE);
329
330   // Register Histograms
331   TH1 * h = NULL;
332   TAxis *ax = NULL;
333   if(!(h = (TH1F *)gROOT->FindObject("hNcls"))){
334     h = new TH1F("hNcls", "N_{clusters} / track", 181, -0.5, 180.5);
335     h->GetXaxis()->SetTitle("N_{clusters}");
336     h->GetYaxis()->SetTitle("Entries");
337   } else h->Reset();
338   fContainer->AddAt(h, kNclustersTrack);
339
340   if(!(h = (TH1F *)gROOT->FindObject("hNclTls"))){
341     h = new TH1F("hNclTls","N_{clusters} / tracklet", 51, -0.5, 50.5);
342     h->GetXaxis()->SetTitle("N_{clusters}");
343     h->GetYaxis()->SetTitle("Entries");
344   } else h->Reset();
345   fContainer->AddAt(h, kNclustersTracklet);
346
347   if(!(h = (TH1F *)gROOT->FindObject("hNtls"))){
348     h = new TH1F("hNtls", "N_{tracklets} / track", AliTRDgeometry::kNlayer, 0.5, 6.5);
349     h->GetXaxis()->SetTitle("N^{tracklet}");
350     h->GetYaxis()->SetTitle("freq. [%]");
351   } else h->Reset();
352   fContainer->AddAt(h, kNtrackletsTrack);
353
354   if(!(h = (TH1F *)gROOT->FindObject("htlsSTA"))){
355     h = new TH1F("hNtlsSTA", "N_{tracklets} / track (Stand Alone)", AliTRDgeometry::kNlayer, 0.5, 6.5);
356     h->GetXaxis()->SetTitle("N^{tracklet}");
357     h->GetYaxis()->SetTitle("freq. [%]");
358   }
359   fContainer->AddAt(h, kNtrackletsSTA);
360
361   if(!(h = (TH1F *)gROOT->FindObject("htlsBAR"))){
362     h = new TH1F("hNtlsBAR", "N_{tracklets} / track (Barrel)", AliTRDgeometry::kNlayer, 0.5, 6.5);
363     h->GetXaxis()->SetTitle("N^{tracklet}");
364     h->GetYaxis()->SetTitle("freq. [%]");
365   }
366   fContainer->AddAt(h, kNtrackletsBAR);
367
368   // 
369   if(!(h = (TH1F *)gROOT->FindObject("hNtlsCross"))){
370     h = new TH1F("hNtlsCross", "N_{tracklets}^{cross} / track", 7, -0.5, 6.5);
371     h->GetXaxis()->SetTitle("n_{row cross}");
372     h->GetYaxis()->SetTitle("freq. [%]");
373   } else h->Reset();
374   fContainer->AddAt(h, kNtrackletsCross);
375
376   if(!(h = (TH1F *)gROOT->FindObject("hNtlsFindable"))){
377     h = new TH1F("hNtlsFindable", "Found/Findable Tracklets" , 101, -0.005, 1.005);
378     h->GetXaxis()->SetTitle("r [a.u]");
379     h->GetYaxis()->SetTitle("Entries");
380   } else h->Reset();
381   fContainer->AddAt(h, kNtrackletsFindable);
382
383   if(!(h = (TH1F *)gROOT->FindObject("hNtrks"))){
384     h = new TH1F("hNtrks", "N_{tracks} / event", 100, 0, 100);
385     h->GetXaxis()->SetTitle("N_{tracks}");
386     h->GetYaxis()->SetTitle("Entries");
387   } else h->Reset();
388   fContainer->AddAt(h, kNtracksEvent);
389
390   if(!(h = (TH1F *)gROOT->FindObject("hNtrksSector"))){
391     h = new TH1F("hNtrksSector", "N_{tracks} / sector", AliTRDgeometry::kNsector, -0.5, 17.5);
392     h->GetXaxis()->SetTitle("sector");
393     h->GetYaxis()->SetTitle("freq. [%]");
394   } else h->Reset();
395   fContainer->AddAt(h, kNtracksSector);
396
397   if(!(h = (TH1F*)gROOT->FindObject("hTrackStatus"))){
398     const Int_t nerr = 7;
399     h = new TH1F("hTrackStatus", "Track Status", nerr, -0.5, nerr-0.5);
400     Char_t *label[nerr] = {"OK", "PROL", "PROP", "AJST", "SNP", "TINI", "UPDT"};
401     ax = h->GetXaxis();
402     for(Int_t ierr=nerr; ierr--;) ax->SetBinLabel(ierr+1, label[ierr]);
403     h->SetYTitle("Relative Error to Good [%]");
404   }
405   fContainer->AddAt(h, kTrackStatus);
406
407   TObjArray *arr = new TObjArray(AliTRDgeometry::kNlayer);
408   arr->SetOwner(kTRUE);  arr->SetName("TrackletStatus");
409   fContainer->AddAt(arr, kTrackletStatus);
410   for(Int_t ily=AliTRDgeometry::kNlayer; ily--;){
411     if(!(h = (TH1F *)gROOT->FindObject(Form("hTrackletStatus%d", ily)))){
412       const Int_t nerr = 8;
413       h = new TH1F(Form("hTrackletStatus%d", ily), "Tracklet status", nerr, -0.5, nerr-0.5);
414       h->SetLineColor(ily+1);
415       Char_t *label[nerr] = {"OK", "Geom", "Bound", "NoCl", "NoAttach", "NoClTr", "NoFit", "Chi2"};
416       ax = h->GetXaxis();
417       for(Int_t ierr=nerr; ierr--;) ax->SetBinLabel(ierr+1, label[ierr]);
418       h->SetYTitle("Relative Error to Good [%]");
419     } else h->Reset();
420     arr->AddAt(h, ily);
421   }
422
423   // <PH> histos
424   arr = new TObjArray(2);
425   arr->SetOwner(kTRUE);  arr->SetName("<PH>");
426   fContainer->AddAt(arr, kPH);
427   if(!(h = (TH1F *)gROOT->FindObject("hPHt"))){
428     h = new TProfile("hPHt", "<PH>", 31, -0.5, 30.5);
429     h->GetXaxis()->SetTitle("Time / 100ns");
430     h->GetYaxis()->SetTitle("<PH> [a.u]");
431   } else h->Reset();
432   arr->AddAt(h, 0);
433   if(!(h = (TH1F *)gROOT->FindObject("hPHx")))
434     h = new TProfile("hPHx", "<PH>", 31, -0.08, 4.88);
435   else h->Reset();
436   arr->AddAt(h, 1);
437
438   // Chi2 histos
439   if(!(h = (TH2S*)gROOT->FindObject("hChi2"))){
440     h = new TH2S("hChi2", "#chi^{2} per track", AliTRDgeometry::kNlayer, .5, AliTRDgeometry::kNlayer+.5, 100, 0, 50);
441     h->SetXTitle("ndf");
442     h->SetYTitle("#chi^{2}/ndf");
443     h->SetZTitle("entries");
444   } else h->Reset();
445   fContainer->AddAt(h, kChi2);
446
447   if(!(h = (TH1F *)gROOT->FindObject("hQcl"))){
448     h = new TH1F("hQcl", "Q_{cluster}", 200, 0, 1200);
449     h->GetXaxis()->SetTitle("Q_{cluster} [a.u.]");
450     h->GetYaxis()->SetTitle("Entries");
451   }else h->Reset();
452   fContainer->AddAt(h, kChargeCluster);
453
454   if(!(h = (TH1F *)gROOT->FindObject("hQtrklt"))){
455     h = new TH1F("hQtrklt", "Q_{tracklet}", 6000, 0, 6000);
456     h->GetXaxis()->SetTitle("Q_{tracklet} [a.u.]");
457     h->GetYaxis()->SetTitle("Entries");
458   }else h->Reset();
459   fContainer->AddAt(h, kChargeTracklet);
460
461
462   if(!(h = (TH1F *)gROOT->FindObject("hEventsTrigger")))
463     h = new TH1F("hEventsTrigger", "Trigger Class", 100, 0, 100);
464   else h->Reset();
465   fContainer->AddAt(h, kNeventsTrigger);
466
467   if(!(h = (TH1F *)gROOT->FindObject("hEventsTriggerTracks")))
468     h = new TH1F("hEventsTriggerTracks", "Trigger Class (Tracks)", 100, 0, 100);
469   else h->Reset();
470   fContainer->AddAt(h, kNeventsTriggerTracks);
471
472   if(!(h = (TH1F *)gROOT->FindObject("hTriggerPurity"))){
473     h = new TH1F("hTriggerPurity", "Trigger Purity", 10, -0.5, 9.5);
474     h->GetXaxis()->SetTitle("Trigger Cluster");
475     h->GetYaxis()->SetTitle("freq.");
476   } else h->Reset();
477   fContainer->AddAt(h, kTriggerPurity);
478
479   return fContainer;
480 }
481
482 /*
483 * Plotting Functions
484 */
485
486 //_______________________________________________________
487 TH1 *AliTRDcheckDET::PlotTrackStatus(const AliTRDtrackV1 *track)
488 {
489 //
490 // Plot the track propagation status. The following errors are defined (see AliTRDtrackV1::ETRDtrackError)
491 //   PROL - track prolongation failure
492 //   PROP - track propagation failure
493 //   AJST - crossing sectors failure
494 //   SNP  - too large bending
495 //   TINI - tracklet initialization failure
496 //   UPDT - track position/covariance update failure 
497 //
498 // Performance plot looks as below:
499 //Begin_Html
500 //<img src="TRD/trackStatus.gif">
501 //End_Html
502 //
503   if(track) fkTrack = track;
504   if(!fkTrack){
505     AliWarning("No Track defined.");
506     return 0x0;
507   }
508   TH1 *h = 0x0;
509   if(!(h = dynamic_cast<TH1F *>(fContainer->At(kTrackStatus)))){
510     AliWarning("No Histogram defined.");
511     return 0x0;
512   }
513   h->Fill(fkTrack->GetStatusTRD());
514   return h;
515 }
516
517 //_______________________________________________________
518 TH1 *AliTRDcheckDET::PlotTrackletStatus(const AliTRDtrackV1 *track)
519 {
520 //
521 // Plot the tracklet propagation status. The following errors are defined for tracklet (see AliTRDtrackV1::ETRDlayerError)
522 //   Geom   - 
523 //   Bound  - tracklet too close to chamber walls
524 //   NoCl   - no clusters in the track roads
525 //   NoAttach - fail to attach clusters
526 //   NoClTr - fail to use clusters for fit
527 //   NoFit  - tracklet fit failled
528 //   Chi2   - chi2 tracklet-track over threshold
529 //
530 // Performance plot looks as below:
531 //Begin_Html
532 //<img src="TRD/trackletStatus.gif">
533 //End_Html
534 //
535   if(track) fkTrack = track;
536   if(!fkTrack){
537     AliWarning("No Track defined.");
538     return 0x0;
539   }
540   TObjArray *arr =0x0;
541   if(!(arr = dynamic_cast<TObjArray*>(fContainer->At(kTrackletStatus)))){
542     AliWarning("Histograms not defined.");
543     return 0x0;
544   }
545
546   TH1 *h = 0x0;
547   for(Int_t ily=AliTRDgeometry::kNlayer; ily--;){
548     if(!(h = dynamic_cast<TH1F*>(arr->At(ily)))){
549       AliWarning(Form("Missing histo for layer %d.", ily));
550       continue;
551     }
552     h->Fill(fkTrack->GetStatusTRD(ily));
553   }
554   return h;
555 }
556
557 //_______________________________________________________
558 TH1 *AliTRDcheckDET::PlotNClustersTracklet(const AliTRDtrackV1 *track){
559   //
560   // Plot the mean number of clusters per tracklet
561   //
562   if(track) fkTrack = track;
563   if(!fkTrack){
564     AliWarning("No Track defined.");
565     return 0x0;
566   }
567   TH1 *h = 0x0;
568   if(!(h = dynamic_cast<TH1F *>(fContainer->At(kNclustersTracklet)))){
569     AliWarning("No Histogram defined.");
570     return 0x0;
571   }
572   AliTRDseedV1 *tracklet = 0x0;
573   for(Int_t itl = 0; itl < AliTRDgeometry::kNlayer; itl++){
574     if(!(tracklet = fkTrack->GetTracklet(itl)) || !tracklet->IsOK()) continue;
575     h->Fill(tracklet->GetN2());
576   }
577   return h;
578 }
579
580 //_______________________________________________________
581 TH1 *AliTRDcheckDET::PlotNClustersTrack(const AliTRDtrackV1 *track){
582   //
583   // Plot the number of clusters in one track
584   //
585   if(track) fkTrack = track;
586   if(!fkTrack){
587     AliWarning("No Track defined.");
588     return 0x0;
589   }
590   TH1 *h = 0x0;
591   if(!(h = dynamic_cast<TH1F *>(fContainer->At(kNclustersTrack)))){
592     AliWarning("No Histogram defined.");
593     return 0x0;
594   }
595   
596   Int_t nclusters = 0;
597   AliTRDseedV1 *tracklet = 0x0;
598   for(Int_t itl = 0; itl < AliTRDgeometry::kNlayer; itl++){
599     if(!(tracklet = fkTrack->GetTracklet(itl)) || !tracklet->IsOK()) continue;
600     nclusters += tracklet->GetN();
601     if(fDebugLevel > 2){
602       Int_t crossing = Int_t(tracklet->IsRowCross());
603       Int_t detector = tracklet->GetDetector();
604       Float_t theta = TMath::ATan(tracklet->GetZref(1));
605       Float_t phi = TMath::ATan(tracklet->GetYref(1));
606       Float_t momentum = 0.;
607       Int_t pdg = 0;
608       Int_t kinkIndex = fkESD ? fkESD->GetKinkIndex() : 0;
609       UShort_t nclsTPC = fkESD ? fkESD->GetTPCncls() : 0;
610       if(fkMC){
611         if(fkMC->GetTrackRef()) momentum = fkMC->GetTrackRef()->P();
612         pdg = fkMC->GetPDG();
613       }
614       (*fDebugStream) << "NClustersTrack"
615         << "Detector="  << detector
616         << "crossing="  << crossing
617         << "momentum="  << momentum
618         << "pdg="                               << pdg
619         << "theta="                     << theta
620         << "phi="                               << phi
621         << "kinkIndex=" << kinkIndex
622         << "TPCncls="           << nclsTPC
623         << "nclusters=" << nclusters
624         << "\n";
625     }
626   }
627   h->Fill(nclusters);
628   return h;
629 }
630
631
632 //_______________________________________________________
633 TH1 *AliTRDcheckDET::PlotNTrackletsTrack(const AliTRDtrackV1 *track){
634   //
635   // Plot the number of tracklets
636   //
637   if(track) fkTrack = track;
638   if(!fkTrack){
639     AliWarning("No Track defined.");
640     return 0x0;
641   }
642   TH1 *h = 0x0, *hMethod = 0x0;
643   if(!(h = dynamic_cast<TH1F *>(fContainer->At(kNtrackletsTrack)))){
644     AliWarning("No Histogram defined.");
645     return 0x0;
646   }
647   Int_t status = fkESD->GetStatus();
648 /*  printf("in/out/refit/pid: TRD[%d|%d|%d|%d]\n", status &AliESDtrack::kTRDin ? 1 : 0, status &AliESDtrack::kTRDout ? 1 : 0, status &AliESDtrack::kTRDrefit ? 1 : 0, status &AliESDtrack::kTRDpid ? 1 : 0);*/
649   if((status & AliESDtrack::kTRDin) != 0){
650     // Full BarrelTrack
651     if(!(hMethod = dynamic_cast<TH1F *>(fContainer->At(kNtrackletsBAR))))
652       AliWarning("Method: Barrel.  Histogram not processed!");
653   } else {
654     // Stand alone Track
655     if(!(hMethod = dynamic_cast<TH1F *>(fContainer->At(kNtrackletsSTA))))
656       AliWarning("Method: StandAlone.  Histogram not processed!");
657   }
658   Int_t nTracklets = fkTrack->GetNumberOfTracklets();
659   h->Fill(nTracklets);
660   hMethod->Fill(nTracklets);
661   if(fDebugLevel > 3){
662     if(nTracklets == 1){
663       // If we have one Tracklet, check in which layer this happens
664       Int_t layer = -1;
665       AliTRDseedV1 *tracklet = 0x0;
666       for(Int_t il = 0; il < AliTRDgeometry::kNlayer; il++){
667         if((tracklet = fkTrack->GetTracklet(il)) && tracklet->IsOK()){layer =  il; break;}
668       }
669       (*fDebugStream) << "NTrackletsTrack"
670         << "Layer=" << layer
671         << "\n";
672     }
673   }
674   return h;
675 }
676
677
678 //_______________________________________________________
679 TH1 *AliTRDcheckDET::PlotNTrackletsRowCross(const AliTRDtrackV1 *track){
680   //
681   // Plot the number of tracklets
682   //
683   if(track) fkTrack = track;
684   if(!fkTrack){
685     AliWarning("No Track defined.");
686     return 0x0;
687   }
688   TH1 *h = 0x0;
689   if(!(h = dynamic_cast<TH1F *>(fContainer->At(kNtrackletsCross)))){
690     AliWarning("No Histogram defined.");
691     return 0x0;
692   }
693
694   Int_t ncross = 0;
695   AliTRDseedV1 *tracklet = 0x0;
696   for(Int_t il = 0; il < AliTRDgeometry::kNlayer; il++){
697     if(!(tracklet = fkTrack->GetTracklet(il)) || !tracklet->IsOK()) continue;
698
699     if(tracklet->IsRowCross()) ncross++;
700   }
701   h->Fill(ncross);
702   return h;
703 }
704
705 //_______________________________________________________
706 TH1 *AliTRDcheckDET::PlotFindableTracklets(const AliTRDtrackV1 *track){
707   //
708   // Plots the ratio of number of tracklets vs.
709   // number of findable tracklets
710   //
711   // Findable tracklets are defined as track prolongation
712   // to layer i does not hit the dead area +- epsilon
713   //
714   // In order to check whether tracklet hist active area in Layer i, 
715   // the track is refitted and the fitted position + an uncertainty 
716   // range is compared to the chamber border (also with a different
717   // uncertainty)
718   //
719   // For the track fit two cases are distinguished:
720   // If the track is a stand alone track (defined by the status bit 
721   // encoding, then the track is fitted with the tilted Rieman model
722   // Otherwise the track is fitted with the Kalman fitter in two steps:
723   // Since the track parameters are give at the outer point, we first 
724   // fit in direction inwards. Afterwards we fit again in direction outwards
725   // to extrapolate the track to layers which are not reached by the track
726   // For the Kalman model, the radial track points have to be shifted by
727   // a distance epsilon in the direction that we want to fit
728   //
729   const Float_t epsilon = 0.01;   // dead area tolerance
730   const Float_t epsilonR = 1;    // shift in radial direction of the anode wire position (Kalman filter only)
731   const Float_t deltaY = 0.7;    // Tolerance in the track position in y-direction
732   const Float_t deltaZ = 7.0;    // Tolerance in the track position in z-direction (Padlength)
733   Double_t xAnode[AliTRDgeometry::kNlayer] = {300.2, 312.8, 325.4, 338.0, 350.6, 363.2}; // Take the default X0
734  
735   if(track) fkTrack = track;
736   if(!fkTrack){
737     AliWarning("No Track defined.");
738     return 0x0;
739   }
740   TH1 *h = 0x0;
741   if(!(h = dynamic_cast<TH1F *>(fContainer->At(kNtrackletsFindable)))){
742     AliWarning("No Histogram defined.");
743     return 0x0;
744   }
745   Int_t nFound = 0, nFindable = 0;
746   Int_t stack = -1;
747   Double_t ymin = 0., ymax = 0., zmin = 0., zmax = 0.;
748   Double_t y = 0., z = 0.;
749   AliTRDseedV1 *tracklet = 0x0;
750   AliTRDpadPlane *pp;  
751   for(Int_t il = 0; il < AliTRDgeometry::kNlayer; il++){
752     if((tracklet = fkTrack->GetTracklet(il)) && tracklet->IsOK()){
753       tracklet->SetReconstructor(fReconstructor);
754       nFound++;
755     }
756   }
757   // 2 Different cases:
758   // 1st stand alone: here we cannot propagate, but be can do a Tilted Rieman Fit
759   // 2nd barrel track: here we propagate the track to the layers
760   AliTrackPoint points[6];
761   Float_t xyz[3];
762   memset(xyz, 0, sizeof(Float_t) * 3);
763   if(((fkESD->GetStatus() & AliESDtrack::kTRDout) > 0) && !((fkESD->GetStatus() & AliESDtrack::kTRDin) > 0)){
764     // stand alone track
765     for(Int_t il = 0; il < AliTRDgeometry::kNlayer; il++){
766       xyz[0] = xAnode[il];
767       points[il].SetXYZ(xyz);
768     }
769     AliTRDtrackerV1::FitRiemanTilt(const_cast<AliTRDtrackV1 *>(fkTrack), 0x0, kTRUE, 6, points);
770   } else {
771     // barrel track
772     //
773     // 2 Steps:
774     // -> Kalman inwards
775     // -> Kalman outwards
776     AliTRDtrackV1 copyTrack(*fkTrack);  // Do Kalman on a (non-constant) copy of the track
777     AliTrackPoint pointsInward[6], pointsOutward[6];
778     for(Int_t il = AliTRDgeometry::kNlayer; il--;){
779       // In order to avoid complications in the Kalman filter if the track points have the same radial
780       // position like the tracklets, we have to shift the radial postion of the anode wire by epsilon
781       // in the direction we want to go
782       // The track points have to be in reverse order for the Kalman Filter inwards
783       xyz[0] = xAnode[AliTRDgeometry::kNlayer - il - 1] - epsilonR;
784       pointsInward[il].SetXYZ(xyz);
785       xyz[0] = xAnode[il] + epsilonR;
786       pointsOutward[il].SetXYZ(xyz);
787     }
788     /*for(Int_t ipt = 0; ipt < AliTRDgeometry::kNlayer; ipt++)
789       printf("%d. X = %f\n", ipt, points[ipt].GetX());*/
790     // Kalman inwards
791     AliTRDtrackerV1::FitKalman(&copyTrack, 0x0, kFALSE, 6, pointsInward);
792     memcpy(points, pointsInward, sizeof(AliTrackPoint) * 6); // Preliminary store the inward results in the Array points
793     // Kalman outwards
794     AliTRDtrackerV1::FitKalman(&copyTrack, 0x0, kTRUE, 6, pointsInward);
795     memcpy(points, pointsOutward, sizeof(AliTrackPoint) * AliTRDgeometry::kNlayer);
796   }
797   for(Int_t il = 0; il < AliTRDgeometry::kNlayer; il++){
798     y = points[il].GetY();
799     z = points[il].GetZ();
800     if((stack = fGeo->GetStack(z, il)) < 0) continue; // Not findable
801     pp = fGeo->GetPadPlane(il, stack);
802     ymin = pp->GetCol0() + epsilon;
803     ymax = pp->GetColEnd() - epsilon; 
804     zmin = pp->GetRowEnd() + epsilon; 
805     zmax = pp->GetRow0() - epsilon;
806     // ignore y-crossing (material)
807     if((z + deltaZ > zmin && z - deltaZ < zmax) && (y + deltaY > ymin && y - deltaY < ymax)) nFindable++;
808       if(fDebugLevel > 3){
809         Double_t posTracklet[2] = {tracklet ? tracklet->GetYfit(0) : 0, tracklet ? tracklet->GetZfit(0) : 0};
810         Int_t hasTracklet = tracklet ? 1 : 0;
811         (*fDebugStream)   << "FindableTracklets"
812           << "layer="     << il
813           << "ytracklet=" << posTracklet[0]
814           << "ytrack="    << y
815           << "ztracklet=" << posTracklet[1]
816           << "ztrack="    << z
817           << "tracklet="  << hasTracklet
818           << "\n";
819       }
820   }
821   
822   h->Fill(nFindable > 0 ? TMath::Min(nFound/static_cast<Double_t>(nFindable), 1.) : 1);
823   if(fDebugLevel > 2) AliInfo(Form("Findable[Found]: %d[%d|%f]", nFindable, nFound, nFound/static_cast<Float_t>(nFindable > 0 ? nFindable : 1)));
824   return h;
825 }
826
827
828 //_______________________________________________________
829 TH1 *AliTRDcheckDET::PlotChi2(const AliTRDtrackV1 *track){
830   //
831   // Plot the chi2 of the track
832   //
833   if(track) fkTrack = track;
834   if(!fkTrack){
835     AliWarning("No Track defined.");
836     return 0x0;
837   }
838   TH1 *h = 0x0;
839   if(!(h = dynamic_cast<TH2S*>(fContainer->At(kChi2)))) {
840     AliWarning("No Histogram defined.");
841     return 0x0;
842   }
843   Int_t n = fkTrack->GetNumberOfTracklets();
844   if(!n) return 0x0;
845
846   h->Fill(n, fkTrack->GetChi2()/n);
847   return h;
848 }
849
850
851 //_______________________________________________________
852 TH1 *AliTRDcheckDET::PlotPHt(const AliTRDtrackV1 *track){
853   //
854   // Plot the average pulse height
855   //
856   if(track) fkTrack = track;
857   if(!fkTrack){
858     AliWarning("No Track defined.");
859     return 0x0;
860   }
861   TProfile *h = 0x0;
862   if(!(h = dynamic_cast<TProfile *>(((TObjArray*)(fContainer->At(kPH)))->At(0)))){
863     AliWarning("No Histogram defined.");
864     return 0x0;
865   }
866   AliTRDseedV1 *tracklet = 0x0;
867   AliTRDcluster *c = 0x0;
868   for(Int_t itl = 0; itl < AliTRDgeometry::kNlayer; itl++){
869     if(!(tracklet = fkTrack->GetTracklet(itl)) || !tracklet->IsOK())continue;
870     Int_t crossing = Int_t(tracklet->IsRowCross());
871     Int_t detector = tracklet->GetDetector();
872     tracklet->ResetClusterIter();
873     while((c = tracklet->NextCluster())){
874       if(!c->IsInChamber()) continue;
875       Int_t localtime        = c->GetLocalTimeBin();
876       Double_t absoluteCharge = TMath::Abs(c->GetQ());
877       h->Fill(localtime, absoluteCharge);
878       if(fDebugLevel > 3){
879         Double_t distance[2];
880         GetDistanceToTracklet(distance, tracklet, c);
881         Float_t theta = TMath::ATan(tracklet->GetZref(1));
882         Float_t phi = TMath::ATan(tracklet->GetYref(1));
883         Float_t momentum = 0.;
884         Int_t pdg = 0;
885         Int_t kinkIndex = fkESD ? fkESD->GetKinkIndex() : 0;
886         UShort_t TPCncls = fkESD ? fkESD->GetTPCncls() : 0;
887         if(fkMC){
888           if(fkMC->GetTrackRef()) momentum = fkMC->GetTrackRef()->P();
889           pdg = fkMC->GetPDG();
890         }
891         (*fDebugStream) << "PHt"
892           << "Detector="        << detector
893           << "crossing="        << crossing
894           << "Timebin="         << localtime
895           << "Charge="          << absoluteCharge
896           << "momentum="        << momentum
897           << "pdg="                             << pdg
898           << "theta="                   << theta
899           << "phi="                             << phi
900           << "kinkIndex="       << kinkIndex
901           << "TPCncls="         << TPCncls
902           << "dy="        << distance[0]
903           << "dz="        << distance[1]
904           << "c.="        << c
905           << "\n";
906       }
907     }
908   }
909   return h;
910 }
911
912 //_______________________________________________________
913 TH1 *AliTRDcheckDET::PlotPHx(const AliTRDtrackV1 *track){
914   //
915   // Plots the average pulse height vs the distance from the anode wire
916   // (plus const anode wire offset)
917   //
918   if(track) fkTrack = track;
919   if(!fkTrack){
920     AliWarning("No Track defined.");
921     return 0x0;
922   }
923   TProfile *h = 0x0;
924   if(!(h = dynamic_cast<TProfile *>(((TObjArray*)(fContainer->At(kPH)))->At(1)))){
925     AliWarning("No Histogram defined.");
926     return 0x0;
927   }
928   Float_t offset = .5*AliTRDgeometry::CamHght();
929   AliTRDseedV1 *tracklet = 0x0;
930   AliTRDcluster *c = 0x0;
931   Double_t distance = 0;
932   Double_t x, y;
933   for(Int_t itl = 0; itl < AliTRDgeometry::kNlayer; itl++){
934     if(!(tracklet = fkTrack->GetTracklet(itl)) || !(tracklet->IsOK())) continue;
935     tracklet->ResetClusterIter();
936     while((c = tracklet->NextCluster())){
937       if(!c->IsInChamber()) continue;
938       x = c->GetX()-AliTRDcluster::GetXcorr(c->GetLocalTimeBin());
939       y = c->GetY()-AliTRDcluster::GetYcorr(AliTRDgeometry::GetLayer(c->GetDetector()), c->GetCenter());
940
941       distance = tracklet->GetX0() - (c->GetX() + 0.3) + offset;
942       h->Fill(distance, TMath::Abs(c->GetQ()));
943     }
944   }  
945   return h;
946 }
947
948 //_______________________________________________________
949 TH1 *AliTRDcheckDET::PlotChargeCluster(const AliTRDtrackV1 *track){
950   //
951   // Plot the cluster charge
952   //
953   if(track) fkTrack = track;
954   if(!fkTrack){
955     AliWarning("No Track defined.");
956     return 0x0;
957   }
958   TH1 *h = 0x0;
959   if(!(h = dynamic_cast<TH1F *>(fContainer->At(kChargeCluster)))){
960     AliWarning("No Histogram defined.");
961     return 0x0;
962   }
963   AliTRDseedV1 *tracklet = 0x0;
964   AliTRDcluster *c = 0x0;
965   for(Int_t itl = 0; itl < AliTRDgeometry::kNlayer; itl++){
966     if(!(tracklet = fkTrack->GetTracklet(itl)) || !tracklet->IsOK())continue;
967     for(Int_t itime = 0; itime < AliTRDtrackerV1::GetNTimeBins(); itime++){
968       if(!(c = tracklet->GetClusters(itime))) continue;
969       h->Fill(c->GetQ());
970     }
971   }
972   return h;
973 }
974
975 //_______________________________________________________
976 TH1 *AliTRDcheckDET::PlotChargeTracklet(const AliTRDtrackV1 *track){
977   //
978   // Plot the charge deposit per chamber
979   //
980   if(track) fkTrack = track;
981   if(!fkTrack){
982     AliWarning("No Track defined.");
983     return 0x0;
984   }
985   TH1 *h = 0x0;
986   if(!(h = dynamic_cast<TH1F *>(fContainer->At(kChargeTracklet)))){
987     AliWarning("No Histogram defined.");
988     return 0x0;
989   }
990   AliTRDseedV1 *tracklet = 0x0;
991   AliTRDcluster *c = 0x0;
992   Double_t qTot = 0;
993   Int_t nTracklets =fkTrack->GetNumberOfTracklets();
994   for(Int_t itl = 0x0; itl < AliTRDgeometry::kNlayer; itl++){
995     if(!(tracklet = fkTrack->GetTracklet(itl)) || !tracklet->IsOK()) continue;
996     qTot = 0.;
997     for(Int_t ic = AliTRDseedV1::kNclusters; ic--;){
998       if(!(c = tracklet->GetClusters(ic))) continue;
999       qTot += TMath::Abs(c->GetQ());
1000     }
1001     h->Fill(qTot);
1002     if(fDebugLevel > 3){
1003       Int_t crossing = (Int_t)tracklet->IsRowCross();
1004       Int_t detector = tracklet->GetDetector();
1005       Float_t theta = TMath::ATan(tracklet->GetZfit(1));
1006       Float_t phi = TMath::ATan(tracklet->GetYfit(1));
1007       Float_t momentum = 0.;
1008       Int_t pdg = 0;
1009       Int_t kinkIndex = fkESD ? fkESD->GetKinkIndex() : 0;
1010       UShort_t nclsTPC = fkESD ? fkESD->GetTPCncls() : 0;
1011       if(fkMC){
1012               if(fkMC->GetTrackRef()) momentum = fkMC->GetTrackRef()->P();
1013         pdg = fkMC->GetPDG();
1014       }
1015       (*fDebugStream) << "ChargeTracklet"
1016         << "Detector="  << detector
1017         << "crossing="  << crossing
1018         << "momentum="  << momentum
1019         << "nTracklets="<< nTracklets
1020         << "pdg="                               << pdg
1021         << "theta="                     << theta
1022         << "phi="                               << phi
1023         << "kinkIndex=" << kinkIndex
1024         << "TPCncls="           << nclsTPC
1025         << "QT="        << qTot
1026         << "\n";
1027     }
1028   }
1029   return h;
1030 }
1031
1032 //_______________________________________________________
1033 TH1 *AliTRDcheckDET::PlotNTracksSector(const AliTRDtrackV1 *track){
1034   //
1035   // Plot the number of tracks per Sector
1036   //
1037   if(track) fkTrack = track;
1038   if(!fkTrack){
1039     AliWarning("No Track defined.");
1040     return 0x0;
1041   }
1042   TH1 *h = 0x0;
1043   if(!(h = dynamic_cast<TH1F *>(fContainer->At(kNtracksSector)))){
1044     AliWarning("No Histogram defined.");
1045     return 0x0;
1046   }
1047
1048   // TODO we should compare with
1049   // sector = Int_t(track->GetAlpha() / AliTRDgeometry::GetAlpha());
1050
1051   AliTRDseedV1 *tracklet = 0x0;
1052   Int_t sector = -1;
1053   for(Int_t itl = 0; itl < AliTRDgeometry::kNlayer; itl++){
1054     if(!(tracklet = fkTrack->GetTracklet(itl)) || !tracklet->IsOK()) continue;
1055     sector = static_cast<Int_t>(tracklet->GetDetector()/AliTRDgeometry::kNdets);
1056     break;
1057   }
1058   h->Fill(sector);
1059   return h;
1060 }
1061
1062
1063 //________________________________________________________
1064 void AliTRDcheckDET::SetRecoParam(AliTRDrecoParam *r)
1065 {
1066
1067   fReconstructor->SetRecoParam(r);
1068 }
1069
1070 //________________________________________________________
1071 void AliTRDcheckDET::GetDistanceToTracklet(Double_t *dist, AliTRDseedV1 * const tracklet, AliTRDcluster * const c)
1072 {
1073   Float_t x = c->GetX();
1074   dist[0] = c->GetY() - tracklet->GetYat(x);
1075   dist[1] = c->GetZ() - tracklet->GetZat(x);
1076 }
1077
1078
1079 //_______________________________________________________
1080 TH1* AliTRDcheckDET::MakePlotChi2()
1081 {
1082 // Plot chi2/track normalized to number of degree of freedom 
1083 // (tracklets) and compare with the theoretical distribution.
1084 // 
1085 // Alex Bercuci <A.Bercuci@gsi.de>
1086
1087   TH2S *h2 = (TH2S*)fContainer->At(kChi2);
1088   TF1 f("fChi2", "[0]*pow(x, [1]-1)*exp(-0.5*x)", 0., 50.);
1089   TLegend *leg = new TLegend(.7,.7,.95,.95);
1090   leg->SetBorderSize(1); leg->SetHeader("Tracklets per Track");
1091   TH1D *h1 = 0x0;
1092   Bool_t kFIRST = kTRUE;
1093   for(Int_t il=1; il<=h2->GetNbinsX(); il++){
1094     h1 = h2->ProjectionY(Form("pyChi2%d", il), il, il);
1095     if(h1->Integral()<50) continue;
1096     h1->Scale(1./h1->Integral());
1097     h1->SetMarkerStyle(7);h1->SetMarkerColor(il);
1098     h1->SetLineColor(il);h1->SetLineStyle(2);
1099     f.SetParameter(1, .5*il);f.SetLineColor(il);
1100     h1->Fit(&f, "QW+", kFIRST ? "pc": "pcsame");
1101     leg->AddEntry(h1, Form("%d", il), "l");
1102     if(kFIRST){
1103       h1->GetXaxis()->SetRangeUser(0., 25.);
1104     }
1105     kFIRST = kFALSE;
1106   }
1107   leg->Draw();
1108   gPad->SetLogy();
1109   return h1;
1110 }
1111
1112
1113 //________________________________________________________
1114 TH1* AliTRDcheckDET::MakePlotNTracklets(){
1115   //
1116   // Make nice bar plot of the number of tracklets in each method
1117   //
1118   TH1F *hBAR = (TH1F *)fContainer->FindObject("hNtlsBAR");
1119   TH1F *hSTA = (TH1F *)fContainer->FindObject("hNtlsSTA");
1120   TH1F *hCON = (TH1F *)fContainer->FindObject("hNtls");
1121   TLegend *leg = new TLegend(0.13, 0.75, 0.39, 0.89);
1122   leg->SetBorderSize(1);
1123   leg->SetFillColor(0);
1124
1125   Float_t scale = hCON->Integral();
1126   hCON->Scale(100./scale);
1127   hCON->SetFillColor(kRed);hCON->SetLineColor(kRed);
1128   hCON->SetBarWidth(0.2);
1129   hCON->SetBarOffset(0.6);
1130   hCON->SetStats(kFALSE);
1131   hCON->GetYaxis()->SetRangeUser(0.,40.);
1132   hCON->GetYaxis()->SetTitleOffset(1.2);
1133   hCON->Draw("bar1"); leg->AddEntry(hCON, "Total", "f");
1134   hCON->SetMaximum(55.);
1135
1136   hBAR->Scale(100./scale);
1137   hBAR->SetFillColor(kGreen);hBAR->SetLineColor(kGreen);
1138   hBAR->SetBarWidth(0.2);
1139   hBAR->SetBarOffset(0.2);
1140   hBAR->SetTitle("");
1141   hBAR->SetStats(kFALSE);
1142   hBAR->GetYaxis()->SetRangeUser(0.,40.);
1143   hBAR->GetYaxis()->SetTitleOffset(1.2);
1144   hBAR->Draw("bar1same"); leg->AddEntry(hBAR, "Barrel", "f");
1145
1146   hSTA->Scale(100./scale);
1147   hSTA->SetFillColor(kBlue);hSTA->SetLineColor(kBlue);
1148   hSTA->SetBarWidth(0.2);
1149   hSTA->SetBarOffset(0.4);
1150   hSTA->SetTitle("");
1151   hSTA->SetStats(kFALSE);
1152   hSTA->GetYaxis()->SetRangeUser(0.,40.);
1153   hSTA->GetYaxis()->SetTitleOffset(1.2);
1154   hSTA->Draw("bar1same"); leg->AddEntry(hSTA, "Stand Alone", "f");
1155   leg->Draw();
1156   gPad->Update();
1157   return hCON;
1158 }
1159
1160 //________________________________________________________
1161 TH1* AliTRDcheckDET::MakePlotPulseHeight(){
1162   //
1163   // Create Plot of the Pluse Height Spectrum
1164   //
1165   TH1 *h, *h1, *h2;
1166   TObjArray *arr = (TObjArray*)fContainer->FindObject("<PH>");
1167   h = (TH1F*)arr->At(0);
1168   h->SetMarkerStyle(24);
1169   h->SetMarkerColor(kBlack);
1170   h->SetLineColor(kBlack);
1171   h->Draw("e1");
1172 //   copy the second histogram in a new one with the same x-dimension as the phs with respect to time
1173   h1 = (TH1F *)arr->At(1);
1174   h2 = new TH1F("hphs1","Average PH", 31, -0.5, 30.5);
1175   for(Int_t ibin = h1->GetXaxis()->GetFirst(); ibin < h1->GetNbinsX(); ibin++) 
1176     h2->SetBinContent(ibin, h1->GetBinContent(ibin));
1177   h2->SetMarkerStyle(22);
1178   h2->SetMarkerColor(kBlue);
1179   h2->SetLineColor(kBlue);
1180   h2->Draw("e1same");
1181   gPad->Update();
1182 //   create axis according to the histogram dimensions of the original second histogram
1183   TGaxis *axis = new TGaxis(gPad->GetUxmin(),
1184                     gPad->GetUymax(),
1185                     gPad->GetUxmax(),
1186                     gPad->GetUymax(),
1187                     -0.08, 4.88, 510,"-L");
1188   axis->SetLineColor(kBlue);
1189   axis->SetLabelColor(kBlue);
1190   axis->SetTextColor(kBlue);
1191   axis->SetTitle("x_{0}-x_{c} [cm]");
1192   axis->Draw();
1193   return h1;
1194 }
1195
1196 //________________________________________________________
1197 Bool_t AliTRDcheckDET::MakeBarPlot(TH1 *histo, Int_t color){
1198   //
1199   // Draw nice bar plots
1200   //
1201   if(!histo->GetEntries()) return kFALSE;
1202   histo->Scale(100./histo->Integral());
1203   histo->SetFillColor(color);
1204   histo->SetBarOffset(.2);
1205   histo->SetBarWidth(.6);
1206   histo->Draw("bar1");
1207   return kTRUE;
1208 }