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