1 ////////////////////////////////////////////////////////////////////////////
4 // Basic checks for tracking and detector performance //
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
13 // - Chi2 distribution for tracks
14 // - Charge distribution per cluster and tracklet
15 // - Number of events and tracks per trigger source
17 // - Track and Tracklet propagation status
20 // Anton Andronic <A.Andronic@gsi.de> //
21 // Alexandru Bercuci <A.Bercuci@gsi.de> //
22 // Markus Fasel <M.Fasel@gsi.de> //
24 ////////////////////////////////////////////////////////////////////////////
35 #include <TGraphErrors.h>
37 #include <TLinearFitter.h>
40 #include <TProfile2D.h>
41 #include <TObjArray.h>
43 #include <TObjString.h>
47 #include <TProfile2D.h>
52 #include "AliTRDcluster.h"
53 #include "AliESDHeader.h"
54 #include "AliESDRun.h"
55 #include "AliESDtrack.h"
56 #include "AliExternalTrackParam.h"
57 #include "AliTRDgeometry.h"
58 #include "AliTRDpadPlane.h"
59 #include "AliTRDSimParam.h"
60 #include "AliTRDseedV1.h"
61 #include "AliTRDtrackV1.h"
62 #include "AliTRDtrackerV1.h"
63 #include "AliTRDReconstructor.h"
64 #include "AliTrackReference.h"
65 #include "AliTrackPointArray.h"
66 #include "AliTracker.h"
67 #include "TTreeStream.h"
69 #include "info/AliTRDtrackInfo.h"
70 #include "info/AliTRDeventInfo.h"
71 #include "AliTRDcheckDET.h"
76 ClassImp(AliTRDcheckDET)
78 //_______________________________________________________
79 AliTRDcheckDET::AliTRDcheckDET():
88 // Default constructor
90 SetNameTitle("TRDcheckDET", "Basic TRD data checker");
93 //_______________________________________________________
94 AliTRDcheckDET::AliTRDcheckDET(char* name):
95 AliTRDrecoTask(name, "Basic TRD data checker")
103 // Default constructor
105 DefineInput(2, AliTRDeventInfo::Class());
107 fReconstructor = new AliTRDReconstructor;
108 fReconstructor->SetRecoParam(AliTRDrecoParam::GetLowFluxParam());
109 fGeo = new AliTRDgeometry;
114 //_______________________________________________________
115 AliTRDcheckDET::~AliTRDcheckDET(){
119 if(fTriggerNames) delete fTriggerNames;
120 delete fReconstructor;
124 //_______________________________________________________
125 void AliTRDcheckDET::UserCreateOutputObjects(){
127 // Create Output Objects
129 AliTRDrecoTask::UserCreateOutputObjects();
130 if(!fTriggerNames) fTriggerNames = new TMap();
133 //_______________________________________________________
134 void AliTRDcheckDET::UserExec(Option_t *opt){
136 // Execution function
137 // Filling TRD quality histos
139 fEventInfo = dynamic_cast<AliTRDeventInfo *>(GetInputData(2));
140 AliDebug(2, Form("EventInfo[%p] Header[%p]", (void*)fEventInfo, (void*)(fEventInfo?fEventInfo->GetEventHeader():NULL)));
142 AliTRDrecoTask::UserExec(opt);
144 Int_t nTracks = 0; // Count the number of tracks per event
145 for(Int_t iti = 0; iti < fTracks->GetEntriesFast(); iti++){
146 if(!fTracks->UncheckedAt(iti)) continue;
147 AliTRDtrackInfo *fTrackInfo = dynamic_cast<AliTRDtrackInfo *>(fTracks->UncheckedAt(iti));
148 if(!fTrackInfo->GetTrack()) continue;
152 dynamic_cast<TH1F *>(fContainer->UncheckedAt(kNtracksEvent))->Fill(nTracks);
154 if(!fEventInfo->GetEventHeader()) return; // For trigger statistics event header is essential
155 Int_t triggermask = fEventInfo->GetEventHeader()->GetTriggerMask();
156 TString triggername = fEventInfo->GetRunInfo()->GetFiredTriggerClasses(triggermask);
157 AliDebug(6, Form("Trigger cluster: %d, Trigger class: %s\n", triggermask, triggername.Data()));
158 dynamic_cast<TH1F *>(fContainer->UncheckedAt(kNeventsTrigger))->Fill(triggermask);
161 dynamic_cast<TH1F *>(fContainer->UncheckedAt(kNeventsTriggerTracks))->Fill(triggermask);
162 if(triggermask <= 20 && !fTriggerNames->FindObject(Form("%d", triggermask))){
163 fTriggerNames->Add(new TObjString(Form("%d", triggermask)), new TObjString(triggername));
164 // also set the label for both histograms
165 TH1 *histo = dynamic_cast<TH1F *>(fContainer->UncheckedAt(kNeventsTriggerTracks));
166 histo->GetXaxis()->SetBinLabel(histo->FindBin(triggermask), triggername);
167 histo = dynamic_cast<TH1F *>(fContainer->UncheckedAt(kNeventsTrigger));
168 histo->GetXaxis()->SetBinLabel(histo->FindBin(triggermask), triggername);
173 //_______________________________________________________
174 Bool_t AliTRDcheckDET::PostProcess(){
176 // Do Postprocessing (for the moment set the number of Reference histograms)
181 // Calculate of the trigger clusters purity
182 h = dynamic_cast<TH1F *>(fContainer->FindObject("hEventsTrigger"));
183 TH1F *h1 = dynamic_cast<TH1F *>(fContainer->FindObject("hEventsTriggerTracks"));
185 Float_t purities[20], val = 0;
186 TString triggernames[20];
187 Int_t nTriggerClasses = 0;
188 for(Int_t ibin = 1; ibin <= h->GetNbinsX(); ibin++){
189 if((val = h1->GetBinContent(ibin))){
190 purities[nTriggerClasses] = val;
191 triggernames[nTriggerClasses] = h1->GetXaxis()->GetBinLabel(ibin);
195 h = dynamic_cast<TH1F *>(fContainer->UncheckedAt(kTriggerPurity));
196 TAxis *ax = h->GetXaxis();
197 for(Int_t itrg = 0; itrg < nTriggerClasses; itrg++){
198 h->Fill(itrg, purities[itrg]);
199 ax->SetBinLabel(itrg+1, triggernames[itrg].Data());
201 ax->SetRangeUser(-0.5, nTriggerClasses+.5);
202 h->GetYaxis()->SetRangeUser(0,1);
205 h=dynamic_cast<TH1F*>(fContainer->At(kTrackStatus));
206 Float_t ok = h->GetBinContent(1);
207 Int_t nerr = h->GetNbinsX();
208 for(Int_t ierr=nerr; ierr--;){
209 h->SetBinContent(ierr+1, ok>0.?1.e2*h->GetBinContent(ierr+1)/ok:0.);
211 h->SetBinContent(1, 0.);
215 TObjArray *arr = dynamic_cast<TObjArray*>(fContainer->UncheckedAt(kTrackletStatus));
216 for(Int_t ily = AliTRDgeometry::kNlayer; ily--;){
217 h=dynamic_cast<TH1F*>(arr->At(ily));
218 Float_t okB = h->GetBinContent(1);
219 Int_t nerrB = h->GetNbinsX();
220 for(Int_t ierr=nerrB; ierr--;){
221 h->SetBinContent(ierr+1, okB>0.?1.e2*h->GetBinContent(ierr+1)/okB:0.);
223 h->SetBinContent(1, 0.);
231 //_______________________________________________________
232 void AliTRDcheckDET::MakeSummary(){
234 // Create summary plots for TRD check DET
235 // This function creates 2 summary plots:
236 // - General Quantities
238 // The function will reuse GetRefFigure
241 TCanvas *cOut = new TCanvas(Form("summary%s1", GetName()), Form("Summary 1 for task %s", GetName()), 1024, 768);
244 // Create figures using GetRefFigure
245 cOut->cd(1); GetRefFigure(kFigNtracksEvent);
246 cOut->cd(2); GetRefFigure(kFigNtracksSector);
247 cOut->cd(3); GetRefFigure(kFigNclustersTrack);
248 cOut->cd(4); GetRefFigure(kFigNclustersTracklet);
249 cOut->cd(5); GetRefFigure(kFigNtrackletsTrack);
250 cOut->cd(6); GetRefFigure(kFigNTrackletsP);
251 cOut->cd(7); GetRefFigure(kFigChargeCluster);
252 cOut->cd(8); GetRefFigure(kFigChargeTracklet);
253 cOut->SaveAs(Form("TRDsummary%s1.gif", GetName()));
257 cOut = new TCanvas(Form("summary%s2", GetName()), Form("Summary 2 for task %s", GetName()), 1024, 512);
258 cOut->cd(); GetRefFigure(kFigPH);
259 cOut->SaveAs(Form("TRDsummary%s2.gif", GetName()));
262 // Third Plot: Mean Number of clusters as function of eta, phi and layer
263 cOut = new TCanvas(Form("summary%s3", GetName()), Form("Summary 3 for task %s", GetName()), 1024, 768);
264 cOut->cd(); MakePlotMeanClustersLayer();
265 cOut->SaveAs(Form("TRDsummary%s3.gif", GetName()));
270 //_______________________________________________________
271 Bool_t AliTRDcheckDET::GetRefFigure(Int_t ifig){
273 // Setting Reference Figures
277 TH1 *h = NULL; TObjArray *arr=NULL;
281 case kFigNclustersTrack:
282 (h=(TH1F*)fContainer->FindObject("hNcls"))->Draw("pl");
283 PutTrendValue("NClustersTrack", h->GetMean());
284 PutTrendValue("NClustersTrackRMS", h->GetRMS());
286 case kFigNclustersTracklet:
287 (h =(TH1F*)fContainer->FindObject("hNclTls"))->Draw("pc");
288 PutTrendValue("NClustersTracklet", h->GetMean());
289 PutTrendValue("NClustersTrackletRMS", h->GetRMS());
291 case kFigNtrackletsTrack:
292 h=MakePlotNTracklets();
294 PutTrendValue("NTrackletsTrack", h->GetMean());
295 PutTrendValue("NTrackletsTrackRMS", h->GetRMS());
298 case kFigNTrackletsP:
299 MakePlotnTrackletsVsP();
301 case kFigNtrackletsCross:
302 h = (TH1F*)fContainer->FindObject("hNtlsCross");
303 if(!MakeBarPlot(h, kRed)) break;
304 PutTrendValue("NTrackletsCross", h->GetMean());
305 PutTrendValue("NTrackletsCrossRMS", h->GetRMS());
307 case kFigNtrackletsFindable:
308 h = (TH1F*)fContainer->FindObject("hNtlsFindable");
309 if(!MakeBarPlot(h, kGreen)) break;
310 PutTrendValue("NTrackletsFindable", h->GetMean());
311 PutTrendValue("NTrackletsFindableRMS", h->GetRMS());
313 case kFigNtracksEvent:
314 (h = (TH1F*)fContainer->FindObject("hNtrks"))->Draw("pl");
315 PutTrendValue("NTracksEvent", h->GetMean());
316 PutTrendValue("NTracksEventRMS", h->GetRMS());
318 case kFigNtracksSector:
319 h = (TH1F*)fContainer->FindObject("hNtrksSector");
320 if(!MakeBarPlot(h, kGreen)) break;
321 PutTrendValue("NTracksSector", h->Integral()/h->GetNbinsX());
323 case kFigTrackStatus:
324 if(!(h=(TH1F *)fContainer->FindObject("hTrackStatus"))) break;
325 h->GetXaxis()->SetRangeUser(0.5, -1);
326 h->GetYaxis()->CenterTitle();
328 PutTrendValue("TrackStatus", h->Integral());
331 case kFigTrackletStatus:
332 if(!(arr = dynamic_cast<TObjArray*>(fContainer->At(kTrackletStatus)))) break;
333 leg = new TLegend(.68, .7, .97, .97);
334 leg->SetBorderSize(0);leg->SetFillStyle(0);
335 leg->SetHeader("TRD layer");
336 for(Int_t ily=AliTRDgeometry::kNlayer; ily--;){
337 if(!(h=dynamic_cast<TH1F*>(arr->At(ily)))) continue;
340 h->GetXaxis()->SetRangeUser(0.5, -1);
341 h->GetYaxis()->CenterTitle();
343 } else h->Draw("samepl");
344 leg->AddEntry(h, Form("ly = %d", ily), "l");
345 PutTrendValue(Form("TrackletStatus%d", ily), h->Integral());
355 gPad->SetMargin(0.125, 0.015, 0.1, 0.1);
356 MakePlotPulseHeight();
359 case kFigChargeCluster:
360 (h = (TH1F*)fContainer->FindObject("hQcl"))->Draw("c");
362 PutTrendValue("ChargeCluster", h->GetMaximumBin());
363 PutTrendValue("ChargeClusterRMS", h->GetRMS());
365 case kFigChargeTracklet:
366 (h=(TH1F*)fContainer->FindObject("hQtrklt"))->Draw("c");
367 PutTrendValue("ChargeTracklet", h->GetMaximumBin());
368 PutTrendValue("ChargeTrackletRMS", h->GetRMS());
370 case kFigNeventsTrigger:
371 ((TH1F*)fContainer->FindObject("hEventsTrigger"))->Draw("");
373 case kFigNeventsTriggerTracks:
374 ((TH1F*)fContainer->FindObject("hEventsTriggerTracks"))->Draw("");
376 case kFigTriggerPurity:
377 if(!MakeBarPlot((TH1F*)fContainer->FindObject("hTriggerPurity"), kGreen)) break;
382 AliInfo(Form("Reference plot [%d] missing result", ifig));
386 //_______________________________________________________
387 TObjArray *AliTRDcheckDET::Histos(){
389 // Create QA histograms
392 if(fContainer) return fContainer;
394 fContainer = new TObjArray(20);
395 //fContainer->SetOwner(kTRUE);
397 // Register Histograms
400 if(!(h = (TH1F *)gROOT->FindObject("hNcls"))){
401 h = new TH1F("hNcls", "N_{clusters} / track", 181, -0.5, 180.5);
402 h->GetXaxis()->SetTitle("N_{clusters}");
403 h->GetYaxis()->SetTitle("Entries");
405 fContainer->AddAt(h, kNclustersTrack);
407 TObjArray *arr = new TObjArray(AliTRDgeometry::kNlayer);
408 arr->SetOwner(kTRUE); arr->SetName("clusters");
409 fContainer->AddAt(arr, kNclustersLayer);
410 for(Int_t ily=AliTRDgeometry::kNlayer; ily--;){
411 if(!(h = (TProfile2D *)gROOT->FindObject(Form("hNcl%d", ily)))){
412 h = new TProfile2D(Form("hNcl%d", ily), Form("Mean Number of clusters in Layer %d", ily), 100, -1.0, 1.0, 50, -1.1*TMath::Pi(), 1.1*TMath::Pi());
413 h->GetXaxis()->SetTitle("#eta");
414 h->GetYaxis()->SetTitle("#phi");
419 if(!(h = (TH1F *)gROOT->FindObject("hNclTls"))){
420 h = new TH1F("hNclTls","N_{clusters} / tracklet", 51, -0.5, 50.5);
421 h->GetXaxis()->SetTitle("N_{clusters}");
422 h->GetYaxis()->SetTitle("Entries");
424 fContainer->AddAt(h, kNclustersTracklet);
426 if(!(h = (TH1F *)gROOT->FindObject("hNtls"))){
427 h = new TH1F("hNtls", "N_{tracklets} / track", AliTRDgeometry::kNlayer, 0.5, 6.5);
428 h->GetXaxis()->SetTitle("N^{tracklet}");
429 h->GetYaxis()->SetTitle("freq. [%]");
431 fContainer->AddAt(h, kNtrackletsTrack);
433 if(!(h = (TH1F *)gROOT->FindObject("htlsSTA"))){
434 h = new TH1F("hNtlsSTA", "N_{tracklets} / track (Stand Alone)", AliTRDgeometry::kNlayer, 0.5, 6.5);
435 h->GetXaxis()->SetTitle("N^{tracklet}");
436 h->GetYaxis()->SetTitle("freq. [%]");
438 fContainer->AddAt(h, kNtrackletsSTA);
440 // Binning for momentum dependent tracklet Plots
443 Float_t binsP[kNp+1], binsTrklt[AliTRDgeometry::kNlayer+1];
444 for(Int_t i=0;i<kNp+1; i++,P+=(TMath::Exp(i*i*.001)-1.)) binsP[i]=P;
445 for(Int_t il = 0; il <= AliTRDgeometry::kNlayer; il++) binsTrklt[il] = 0.5 + il;
446 if(!(h = (TH1F *)gROOT->FindObject("htlsBAR"))){
447 // Make tracklets for barrel tracks momentum dependent (if we do not exceed min and max values)
448 h = new TH2F("hNtlsBAR",
449 "N_{tracklets} / track;p [GeV/c];N^{tracklet};freq. [%]",
450 kNp, binsP, AliTRDgeometry::kNlayer, binsTrklt);
452 fContainer->AddAt(h, kNtrackletsBAR);
455 if(!(h = (TH1F *)gROOT->FindObject("hNtlsCross"))){
456 h = new TH1F("hNtlsCross", "N_{tracklets}^{cross} / track", 7, -0.5, 6.5);
457 h->GetXaxis()->SetTitle("n_{row cross}");
458 h->GetYaxis()->SetTitle("freq. [%]");
460 fContainer->AddAt(h, kNtrackletsCross);
462 if(!(h = (TH1F *)gROOT->FindObject("hNtlsFindable"))){
463 h = new TH1F("hNtlsFindable", "Found/Findable Tracklets" , 101, -0.005, 1.005);
464 h->GetXaxis()->SetTitle("r [a.u]");
465 h->GetYaxis()->SetTitle("Entries");
467 fContainer->AddAt(h, kNtrackletsFindable);
469 if(!(h = (TH1F *)gROOT->FindObject("hNtrks"))){
470 h = new TH1F("hNtrks", "N_{tracks} / event", 100, 0, 100);
471 h->GetXaxis()->SetTitle("N_{tracks}");
472 h->GetYaxis()->SetTitle("Entries");
474 fContainer->AddAt(h, kNtracksEvent);
476 if(!(h = (TH1F *)gROOT->FindObject("hNtrksSector"))){
477 h = new TH1F("hNtrksSector", "N_{tracks} / sector", AliTRDgeometry::kNsector, -0.5, 17.5);
478 h->GetXaxis()->SetTitle("sector");
479 h->GetYaxis()->SetTitle("freq. [%]");
481 fContainer->AddAt(h, kNtracksSector);
483 if(!(h = (TH1F*)gROOT->FindObject("hTrackStatus"))){
484 const Int_t nerr = 7;
485 h = new TH1F("hTrackStatus", "Track Status", nerr, -0.5, nerr-0.5);
486 const Char_t *label[nerr] = {"OK", "PROL", "PROP", "AJST", "SNP", "TINI", "UPDT"};
488 for(Int_t ierr=nerr; ierr--;) ax->SetBinLabel(ierr+1, label[ierr]);
489 h->SetYTitle("Relative Error to Good [%]");
491 fContainer->AddAt(h, kTrackStatus);
493 arr = new TObjArray(AliTRDgeometry::kNlayer);
494 arr->SetOwner(kTRUE); arr->SetName("TrackletStatus");
495 fContainer->AddAt(arr, kTrackletStatus);
496 for(Int_t ily=AliTRDgeometry::kNlayer; ily--;){
497 if(!(h = (TH1F *)gROOT->FindObject(Form("hTrackletStatus%d", ily)))){
498 const Int_t nerr = 8;
499 h = new TH1F(Form("hTrackletStatus%d", ily), "Tracklet status", nerr, -0.5, nerr-0.5);
500 h->SetLineColor(ily+1);
501 const Char_t *label[nerr] = {"OK", "Geom", "Bound", "NoCl", "NoAttach", "NoClTr", "NoFit", "Chi2"};
503 for(Int_t ierr=nerr; ierr--;) ax->SetBinLabel(ierr+1, label[ierr]);
504 h->SetYTitle("Relative Error to Good [%]");
510 arr = new TObjArray(3);
511 arr->SetOwner(kTRUE); arr->SetName("<PH>");
512 fContainer->AddAt(arr, kPH);
513 if(!(h = (TH1F *)gROOT->FindObject("hPHt"))){
514 h = new TProfile("hPHt", "<PH>", 31, -0.5, 30.5);
515 h->GetXaxis()->SetTitle("Time / 100ns");
516 h->GetYaxis()->SetTitle("<PH> [a.u]");
519 if(!(h = (TH1F *)gROOT->FindObject("hPHx")))
520 h = new TProfile("hPHx", "<PH>", 31, -0.08, 4.88);
523 if(!(h = (TH2F *)gROOT->FindObject("hPH2D"))){
524 h = new TH2F("hPH2D", "Charge Distribution / time", 31, -0.5, 30.5, 100, 0, 1024);
525 h->GetXaxis()->SetTitle("Time / 100ns");
526 h->GetYaxis()->SetTitle("Charge / a.u.");
531 if(!(h = (TH2S*)gROOT->FindObject("hChi2"))){
532 h = new TH2S("hChi2", "#chi^{2} per track", AliTRDgeometry::kNlayer, .5, AliTRDgeometry::kNlayer+.5, 100, 0, 50);
534 h->SetYTitle("#chi^{2}/ndf");
535 h->SetZTitle("entries");
537 fContainer->AddAt(h, kChi2);
539 if(!(h = (TH1F *)gROOT->FindObject("hQcl"))){
540 h = new TH1F("hQcl", "Q_{cluster}", 200, 0, 1200);
541 h->GetXaxis()->SetTitle("Q_{cluster} [a.u.]");
542 h->GetYaxis()->SetTitle("Entries");
544 fContainer->AddAt(h, kChargeCluster);
546 if(!(h = (TH1F *)gROOT->FindObject("hQtrklt"))){
547 h = new TH1F("hQtrklt", "Q_{tracklet}", 6000, 0, 6000);
548 h->GetXaxis()->SetTitle("Q_{tracklet} [a.u.]");
549 h->GetYaxis()->SetTitle("Entries");
551 fContainer->AddAt(h, kChargeTracklet);
554 if(!(h = (TH1F *)gROOT->FindObject("hEventsTrigger")))
555 h = new TH1F("hEventsTrigger", "Trigger Class", 100, 0, 100);
557 printf("Histos Adding \n");
559 fContainer->AddAt(h, kNeventsTrigger);
561 if(!(h = (TH1F *)gROOT->FindObject("hEventsTriggerTracks")))
562 h = new TH1F("hEventsTriggerTracks", "Trigger Class (Tracks)", 100, 0, 100);
564 fContainer->AddAt(h, kNeventsTriggerTracks);
566 if(!(h = (TH1F *)gROOT->FindObject("hTriggerPurity"))){
567 h = new TH1F("hTriggerPurity", "Trigger Purity", 10, -0.5, 9.5);
568 h->GetXaxis()->SetTitle("Trigger Cluster");
569 h->GetYaxis()->SetTitle("freq.");
571 fContainer->AddAt(h, kTriggerPurity);
580 //_______________________________________________________
581 TH1 *AliTRDcheckDET::PlotTrackStatus(const AliTRDtrackV1 *track)
584 // Plot the track propagation status. The following errors are defined (see AliTRDtrackV1::ETRDtrackError)
585 // PROL - track prolongation failure
586 // PROP - track propagation failure
587 // AJST - crossing sectors failure
588 // SNP - too large bending
589 // TINI - tracklet initialization failure
590 // UPDT - track position/covariance update failure
592 // Performance plot looks as below:
594 //<img src="TRD/trackStatus.gif">
597 if(track) fkTrack = track;
599 AliDebug(4, "No Track defined.");
603 if(!(h = dynamic_cast<TH1F *>(fContainer->At(kTrackStatus)))){
604 AliWarning("No Histogram defined.");
607 h->Fill(fkTrack->GetStatusTRD());
611 //_______________________________________________________
612 TH1 *AliTRDcheckDET::PlotTrackletStatus(const AliTRDtrackV1 *track)
615 // Plot the tracklet propagation status. The following errors are defined for tracklet (see AliTRDtrackV1::ETRDlayerError)
617 // Bound - tracklet too close to chamber walls
618 // NoCl - no clusters in the track roads
619 // NoAttach - fail to attach clusters
620 // NoClTr - fail to use clusters for fit
621 // NoFit - tracklet fit failled
622 // Chi2 - chi2 tracklet-track over threshold
624 // Performance plot looks as below:
626 //<img src="TRD/trackletStatus.gif">
629 if(track) fkTrack = track;
631 AliDebug(4, "No Track defined.");
634 TObjArray *arr =NULL;
635 if(!(arr = dynamic_cast<TObjArray*>(fContainer->At(kTrackletStatus)))){
636 AliWarning("Histograms not defined.");
641 for(Int_t ily=AliTRDgeometry::kNlayer; ily--;){
642 if(!(h = dynamic_cast<TH1F*>(arr->At(ily)))){
643 AliWarning(Form("Missing histo for layer %d.", ily));
646 h->Fill(fkTrack->GetStatusTRD(ily));
651 //_______________________________________________________
652 TH1 *AliTRDcheckDET::PlotNClustersTracklet(const AliTRDtrackV1 *track){
654 // Plot the mean number of clusters per tracklet
656 if(track) fkTrack = track;
658 AliDebug(4, "No Track defined.");
661 AliExternalTrackParam *par = fkTrack->GetTrackIn() ? fkTrack->GetTrackIn() : fkTrack->GetTrackOut();
663 TProfile2D *hlayer = NULL;
664 Double_t eta = 0., phi = 0.;
665 if(!(h = dynamic_cast<TH1F *>(fContainer->At(kNclustersTracklet)))){
666 AliWarning("No Histogram defined.");
669 AliTRDseedV1 *tracklet = NULL;
670 TObjArray *histosLayer = dynamic_cast<TObjArray *>(fContainer->At(kNclustersLayer));
672 AliWarning("No Histograms for single layer defined");
674 for(Int_t itl = 0; itl < AliTRDgeometry::kNlayer; itl++){
675 if(!(tracklet = fkTrack->GetTracklet(itl)) || !tracklet->IsOK()) continue;
676 h->Fill(tracklet->GetN2());
677 if(histosLayer && par){
678 if((hlayer = dynamic_cast<TProfile2D *>(histosLayer->At(itl)))){
679 GetEtaPhiAt(par, tracklet->GetX0(), eta, phi);
680 hlayer->Fill(eta, phi, tracklet->GetN2());
687 //_______________________________________________________
688 TH1 *AliTRDcheckDET::PlotNClustersTrack(const AliTRDtrackV1 *track){
690 // Plot the number of clusters in one track
692 if(track) fkTrack = track;
694 AliDebug(4, "No Track defined.");
698 if(!(h = dynamic_cast<TH1F *>(fContainer->At(kNclustersTrack)))){
699 AliWarning("No Histogram defined.");
704 AliTRDseedV1 *tracklet = NULL;
705 AliExternalTrackParam *par = fkTrack->GetTrackOut() ? fkTrack->GetTrackOut() : fkTrack->GetTrackIn();
706 if(!par) return NULL;
707 Double_t momentumRec = par->P();
708 for(Int_t itl = 0; itl < AliTRDgeometry::kNlayer; itl++){
709 if(!(tracklet = fkTrack->GetTracklet(itl)) || !tracklet->IsOK()) continue;
710 Int_t n(tracklet->GetN());
712 if(DebugLevel() > 2){
713 Int_t crossing = Int_t(tracklet->IsRowCross());
714 Int_t detector = tracklet->GetDetector();
715 Float_t theta = TMath::ATan(tracklet->GetZref(1));
716 Float_t phi = TMath::ATan(tracklet->GetYref(1));
717 Float_t momentumMC = 0.;
719 Int_t kinkIndex = fkESD ? fkESD->GetKinkIndex() : 0;
720 UShort_t nclsTPC = fkESD ? fkESD->GetTPCncls() : 0;
722 if(fkMC->GetTrackRef()) momentumMC = fkMC->GetTrackRef()->P();
723 pdg = fkMC->GetPDG();
725 (*DebugStream()) << "NClustersTrack"
726 << "Detector=" << detector
727 << "crossing=" << crossing
728 << "momentumMC="<< momentumMC
729 << "momentumRec=" << momentumRec
733 << "kinkIndex=" << kinkIndex
734 << "TPCncls=" << nclsTPC
744 //_______________________________________________________
745 TH1 *AliTRDcheckDET::PlotNTrackletsTrack(const AliTRDtrackV1 *track){
747 // Plot the number of tracklets
749 if(track) fkTrack = track;
751 AliDebug(4, "No Track defined.");
754 TH1 *h = NULL, *hSta = NULL; TH2 *hBarrel = NULL;
755 if(!(h = dynamic_cast<TH1F *>(fContainer->At(kNtrackletsTrack)))){
756 AliWarning("No Histogram defined.");
759 Int_t nTracklets = fkTrack->GetNumberOfTracklets();
762 Int_t status = fkESD->GetStatus();
764 /* 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);*/
766 Int_t method = -1; // to distinguish between stand alone and full barrel tracks in the debugging
767 if((status & AliESDtrack::kTRDin) != 0){
769 // Full Barrel Track: Save momentum dependence
770 if(!(hBarrel = dynamic_cast<TH2F *>(fContainer->At(kNtrackletsBAR)))){
771 AliWarning("Method: Barrel. Histogram not processed!");
774 AliExternalTrackParam *par(fkTrack->GetTrackIn());
776 AliError("Input track params missing");
779 p = par->P(); // p needed later in the debug streaming
780 hBarrel->Fill(p, nTracklets);
782 // Stand alone Track: momentum dependence not usefull
784 if(!(hSta = dynamic_cast<TH1F *>(fContainer->At(kNtrackletsSTA)))) {
785 AliWarning("Method: StandAlone. Histogram not processed!");
788 hSta->Fill(nTracklets);
791 if(DebugLevel() > 2){
792 AliTRDseedV1 *tracklet = NULL;
793 Int_t sector = -1, stack = -1, detector;
794 for(Int_t itl = 0; itl < AliTRDgeometry::kNlayer; itl++){
795 if(!(tracklet = fkTrack->GetTracklet(itl)) || !(tracklet->IsOK())) continue;
796 detector = tracklet->GetDetector();
797 sector = fGeo->GetSector(detector);
798 stack = fGeo->GetStack(detector);
801 (*DebugStream()) << "NTrackletsTrack"
802 << "Sector=" << sector
804 << "NTracklets=" << nTracklets
805 << "Method=" << method
809 if(DebugLevel() > 3){
810 AliTRDseedV1 *tracklet = NULL;
811 for(Int_t il = 0; il < AliTRDgeometry::kNlayer; il++){
812 if((tracklet = fkTrack->GetTracklet(il)) && tracklet->IsOK()){
813 (*DebugStream()) << "NTrackletsLayer"
824 //_______________________________________________________
825 TH1 *AliTRDcheckDET::PlotNTrackletsRowCross(const AliTRDtrackV1 *track){
827 // Plot the number of tracklets
829 if(track) fkTrack = track;
831 AliDebug(4, "No Track defined.");
835 if(!(h = dynamic_cast<TH1F *>(fContainer->At(kNtrackletsCross)))){
836 AliWarning("No Histogram defined.");
841 AliTRDseedV1 *tracklet = NULL;
842 for(Int_t il = 0; il < AliTRDgeometry::kNlayer; il++){
843 if(!(tracklet = fkTrack->GetTracklet(il)) || !tracklet->IsOK()) continue;
845 if(tracklet->IsRowCross()) ncross++;
851 //_______________________________________________________
852 TH1 *AliTRDcheckDET::PlotFindableTracklets(const AliTRDtrackV1 *track){
854 // Plots the ratio of number of tracklets vs.
855 // number of findable tracklets
857 // Findable tracklets are defined as track prolongation
858 // to layer i does not hit the dead area +- epsilon
860 // In order to check whether tracklet hist active area in Layer i,
861 // the track is refitted and the fitted position + an uncertainty
862 // range is compared to the chamber border (also with a different
865 // For the track fit two cases are distinguished:
866 // If the track is a stand alone track (defined by the status bit
867 // encoding, then the track is fitted with the tilted Rieman model
868 // Otherwise the track is fitted with the Kalman fitter in two steps:
869 // Since the track parameters are give at the outer point, we first
870 // fit in direction inwards. Afterwards we fit again in direction outwards
871 // to extrapolate the track to layers which are not reached by the track
872 // For the Kalman model, the radial track points have to be shifted by
873 // a distance epsilon in the direction that we want to fit
875 const Float_t epsilon = 0.01; // dead area tolerance
876 const Float_t epsilonR = 1; // shift in radial direction of the anode wire position (Kalman filter only)
877 const Float_t deltaY = 0.7; // Tolerance in the track position in y-direction
878 const Float_t deltaZ = 7.0; // Tolerance in the track position in z-direction (Padlength)
879 Double_t xAnode[AliTRDgeometry::kNlayer] = {300.2, 312.8, 325.4, 338.0, 350.6, 363.2}; // Take the default X0
881 if(track) fkTrack = track;
883 AliDebug(4, "No Track defined.");
887 if(!(h = dynamic_cast<TH1F *>(fContainer->At(kNtrackletsFindable)))){
888 AliWarning("No Histogram defined.");
891 Int_t nFound = 0, nFindable = 0;
893 Double_t ymin = 0., ymax = 0., zmin = 0., zmax = 0.;
894 Double_t y = 0., z = 0.;
895 AliTRDseedV1 *tracklet = NULL;
897 for(Int_t il = 0; il < AliTRDgeometry::kNlayer; il++){
898 if((tracklet = fkTrack->GetTracklet(il)) && tracklet->IsOK()){
899 tracklet->SetReconstructor(fReconstructor);
903 // 2 Different cases:
904 // 1st stand alone: here we cannot propagate, but be can do a Tilted Rieman Fit
905 // 2nd barrel track: here we propagate the track to the layers
906 AliTrackPoint points[6];
908 memset(xyz, 0, sizeof(Float_t) * 3);
909 if(((fkESD->GetStatus() & AliESDtrack::kTRDout) > 0) && !((fkESD->GetStatus() & AliESDtrack::kTRDin) > 0)){
911 for(Int_t il = 0; il < AliTRDgeometry::kNlayer; il++){
913 points[il].SetXYZ(xyz);
915 AliTRDtrackerV1::FitRiemanTilt(const_cast<AliTRDtrackV1 *>(fkTrack), NULL, kTRUE, 6, points);
921 // -> Kalman outwards
922 AliTRDtrackV1 copyTrack(*fkTrack); // Do Kalman on a (non-constant) copy of the track
923 AliTrackPoint pointsInward[6], pointsOutward[6];
924 for(Int_t il = AliTRDgeometry::kNlayer; il--;){
925 // In order to avoid complications in the Kalman filter if the track points have the same radial
926 // position like the tracklets, we have to shift the radial postion of the anode wire by epsilon
927 // in the direction we want to go
928 // The track points have to be in reverse order for the Kalman Filter inwards
929 xyz[0] = xAnode[AliTRDgeometry::kNlayer - il - 1] - epsilonR;
930 pointsInward[il].SetXYZ(xyz);
931 xyz[0] = xAnode[il] + epsilonR;
932 pointsOutward[il].SetXYZ(xyz);
934 /*for(Int_t ipt = 0; ipt < AliTRDgeometry::kNlayer; ipt++)
935 printf("%d. X = %f\n", ipt, points[ipt].GetX());*/
937 AliTRDtrackerV1::FitKalman(©Track, NULL, kFALSE, 6, pointsInward);
938 memcpy(points, pointsInward, sizeof(AliTrackPoint) * 6); // Preliminary store the inward results in the Array points
940 AliTRDtrackerV1::FitKalman(©Track, NULL, kTRUE, 6, pointsInward);
941 memcpy(points, pointsOutward, sizeof(AliTrackPoint) * AliTRDgeometry::kNlayer);
943 for(Int_t il = 0; il < AliTRDgeometry::kNlayer; il++){
944 y = points[il].GetY();
945 z = points[il].GetZ();
946 if((stack = fGeo->GetStack(z, il)) < 0) continue; // Not findable
947 pp = fGeo->GetPadPlane(il, stack);
948 ymin = pp->GetCol0() + epsilon;
949 ymax = pp->GetColEnd() - epsilon;
950 zmin = pp->GetRowEnd() + epsilon;
951 zmax = pp->GetRow0() - epsilon;
952 // ignore y-crossing (material)
953 if((z + deltaZ > zmin && z - deltaZ < zmax) && (y + deltaY > ymin && y - deltaY < ymax)) nFindable++;
954 if(DebugLevel() > 3){
955 Double_t posTracklet[2] = {tracklet ? tracklet->GetYfit(0) : 0, tracklet ? tracklet->GetZfit(0) : 0};
956 Int_t hasTracklet = tracklet ? 1 : 0;
957 (*DebugStream()) << "FindableTracklets"
959 << "ytracklet=" << posTracklet[0]
961 << "ztracklet=" << posTracklet[1]
963 << "tracklet=" << hasTracklet
968 h->Fill(nFindable > 0 ? TMath::Min(nFound/static_cast<Double_t>(nFindable), 1.) : 1);
969 AliDebug(2, Form("Findable[Found]: %d[%d|%f]", nFindable, nFound, nFound/static_cast<Float_t>(nFindable > 0 ? nFindable : 1)));
974 //_______________________________________________________
975 TH1 *AliTRDcheckDET::PlotChi2(const AliTRDtrackV1 *track){
977 // Plot the chi2 of the track
979 if(track) fkTrack = track;
981 AliDebug(4, "No Track defined.");
985 if(!(h = dynamic_cast<TH2S*>(fContainer->At(kChi2)))) {
986 AliWarning("No Histogram defined.");
989 Int_t n = fkTrack->GetNumberOfTracklets();
992 h->Fill(n, fkTrack->GetChi2()/n);
997 //_______________________________________________________
998 TH1 *AliTRDcheckDET::PlotPHt(const AliTRDtrackV1 *track){
1000 // Plot the average pulse height
1002 if(track) fkTrack = track;
1004 AliDebug(4, "No Track defined.");
1007 TProfile *h = NULL; TH2F *phs2D = NULL;
1008 if(!(h = dynamic_cast<TProfile *>(((TObjArray*)(fContainer->At(kPH)))->At(0)))){
1009 AliWarning("No Histogram defined.");
1012 if(!(phs2D = dynamic_cast<TH2F *>(((TObjArray*)(fContainer->At(kPH)))->At(2)))){
1013 AliWarning("2D Pulse Height histogram not defined. Histogramm cannot be filled");
1015 AliTRDseedV1 *tracklet = NULL;
1016 AliTRDcluster *c = NULL;
1017 for(Int_t itl = 0; itl < AliTRDgeometry::kNlayer; itl++){
1018 if(!(tracklet = fkTrack->GetTracklet(itl)) || !tracklet->IsOK())continue;
1019 Int_t crossing = Int_t(tracklet->IsRowCross());
1020 Int_t detector = tracklet->GetDetector();
1021 tracklet->ResetClusterIter();
1022 while((c = tracklet->NextCluster())){
1023 if(!IsUsingClustersOutsideChamber() && !c->IsInChamber()) continue;
1024 Int_t localtime = c->GetLocalTimeBin();
1025 Double_t absoluteCharge = TMath::Abs(c->GetQ());
1026 h->Fill(localtime, absoluteCharge);
1027 phs2D->Fill(localtime, absoluteCharge);
1028 if(DebugLevel() > 3){
1029 Int_t inChamber = c->IsInChamber() ? 1 : 0;
1030 Double_t distance[2];
1031 GetDistanceToTracklet(distance, tracklet, c);
1032 Float_t theta = TMath::ATan(tracklet->GetZref(1));
1033 Float_t phi = TMath::ATan(tracklet->GetYref(1));
1034 AliExternalTrackParam *trdPar = fkTrack->GetTrackIn();
1035 Float_t momentumMC = 0, momentumRec = trdPar ? trdPar->P() : track->P(); // prefer Track Low
1037 Int_t kinkIndex = fkESD ? fkESD->GetKinkIndex() : 0;
1038 UShort_t TPCncls = fkESD ? fkESD->GetTPCncls() : 0;
1040 if(fkMC->GetTrackRef()) momentumMC = fkMC->GetTrackRef()->P();
1041 pdg = fkMC->GetPDG();
1043 (*DebugStream()) << "PHt"
1044 << "Detector=" << detector
1045 << "crossing=" << crossing
1046 << "inChamber=" << inChamber
1047 << "Timebin=" << localtime
1048 << "Charge=" << absoluteCharge
1049 << "momentumMC=" << momentumMC
1050 << "momentumRec=" << momentumRec
1052 << "theta=" << theta
1054 << "kinkIndex=" << kinkIndex
1055 << "TPCncls=" << TPCncls
1056 << "dy=" << distance[0]
1057 << "dz=" << distance[1]
1066 //_______________________________________________________
1067 TH1 *AliTRDcheckDET::PlotPHx(const AliTRDtrackV1 *track){
1069 // Plots the average pulse height vs the distance from the anode wire
1070 // (plus const anode wire offset)
1072 if(track) fkTrack = track;
1074 AliDebug(4, "No Track defined.");
1078 if(!(h = dynamic_cast<TProfile *>(((TObjArray*)(fContainer->At(kPH)))->At(1)))){
1079 AliWarning("No Histogram defined.");
1082 AliTRDseedV1 *tracklet(NULL);
1083 AliTRDcluster *c(NULL);
1084 Double_t xd(0.), dqdl(0.);
1085 TVectorD vq(AliTRDseedV1::kNtb), vxd(AliTRDseedV1::kNtb), vdqdl(AliTRDseedV1::kNtb);
1086 for(Int_t itl = 0; itl < AliTRDgeometry::kNlayer; itl++){
1087 if(!(tracklet = fkTrack->GetTracklet(itl)) || !(tracklet->IsOK())) continue;
1088 Int_t det(tracklet->GetDetector());
1089 Bool_t rc(tracklet->IsRowCross());
1090 for(Int_t ic(0); ic<AliTRDseedV1::kNtb; ic++){
1091 Bool_t kFIRST(kFALSE);
1092 if(!(c = tracklet->GetClusters(ic))){
1093 if(!(c = tracklet->GetClusters(AliTRDseedV1::kNtb+ic))) continue;
1094 } else kFIRST=kTRUE;
1095 if(!IsUsingClustersOutsideChamber() && !c->IsInChamber()) continue;
1096 xd = tracklet->GetX0() - c->GetX(); vxd[ic] = xd;
1097 dqdl=tracklet->GetdQdl(ic); vdqdl[ic] = dqdl;
1099 if(kFIRST && (c = tracklet->GetClusters(AliTRDseedV1::kNtb+ic))) vq[ic]+=c->GetQ();
1102 if(DebugLevel() > 3){
1103 (*DebugStream()) << "PHx"
1108 << "dqdl=" << &vdqdl
1115 //_______________________________________________________
1116 TH1 *AliTRDcheckDET::PlotChargeCluster(const AliTRDtrackV1 *track){
1118 // Plot the cluster charge
1120 if(track) fkTrack = track;
1122 AliDebug(4, "No Track defined.");
1126 if(!(h = dynamic_cast<TH1F *>(fContainer->At(kChargeCluster)))){
1127 AliWarning("No Histogram defined.");
1130 AliTRDseedV1 *tracklet = NULL;
1131 AliTRDcluster *c = NULL;
1132 for(Int_t itl = 0; itl < AliTRDgeometry::kNlayer; itl++){
1133 if(!(tracklet = fkTrack->GetTracklet(itl)) || !tracklet->IsOK())continue;
1134 for(Int_t ic(0); ic < AliTRDseedV1::kNtb; ic++){
1135 Bool_t kFIRST(kFALSE);
1136 if(!(c = tracklet->GetClusters(ic))) {
1137 if(!(c = tracklet->GetClusters(AliTRDseedV1::kNtb+ic))) continue;
1138 } else kFIRST = kTRUE;
1139 Float_t q(c->GetQ());
1140 if(kFIRST && (c = tracklet->GetClusters(AliTRDseedV1::kNtb+ic))) q+=c->GetQ();
1147 //_______________________________________________________
1148 TH1 *AliTRDcheckDET::PlotChargeTracklet(const AliTRDtrackV1 *track){
1150 // Plot the charge deposit per chamber
1152 if(track) fkTrack = track;
1154 AliDebug(4, "No Track defined.");
1158 if(!(h = dynamic_cast<TH1F *>(fContainer->At(kChargeTracklet)))){
1159 AliWarning("No Histogram defined.");
1162 AliTRDseedV1 *tracklet = NULL;
1163 AliTRDcluster *c = NULL;
1165 Int_t nTracklets =fkTrack->GetNumberOfTracklets();
1166 for(Int_t itl(0); itl < AliTRDgeometry::kNlayer; itl++){
1167 if(!(tracklet = fkTrack->GetTracklet(itl)) || !tracklet->IsOK()) continue;
1169 for(Int_t ic = AliTRDseedV1::kNclusters; ic--;){
1170 if(!(c = tracklet->GetClusters(ic))) continue;
1171 qTot += TMath::Abs(c->GetQ());
1174 if(DebugLevel() > 3){
1175 Int_t crossing = (Int_t)tracklet->IsRowCross();
1176 Int_t detector = tracklet->GetDetector();
1177 Float_t theta = TMath::ATan(tracklet->GetZfit(1));
1178 Float_t phi = TMath::ATan(tracklet->GetYfit(1));
1179 Float_t momentum = 0.;
1181 Int_t kinkIndex = fkESD ? fkESD->GetKinkIndex() : 0;
1182 UShort_t nclsTPC = fkESD ? fkESD->GetTPCncls() : 0;
1184 if(fkMC->GetTrackRef()) momentum = fkMC->GetTrackRef()->P();
1185 pdg = fkMC->GetPDG();
1187 (*DebugStream()) << "ChargeTracklet"
1188 << "Detector=" << detector
1189 << "crossing=" << crossing
1190 << "momentum=" << momentum
1191 << "nTracklets="<< nTracklets
1193 << "theta=" << theta
1195 << "kinkIndex=" << kinkIndex
1196 << "TPCncls=" << nclsTPC
1204 //_______________________________________________________
1205 TH1 *AliTRDcheckDET::PlotNTracksSector(const AliTRDtrackV1 *track){
1207 // Plot the number of tracks per Sector
1209 if(track) fkTrack = track;
1211 AliDebug(4, "No Track defined.");
1215 if(!(h = dynamic_cast<TH1F *>(fContainer->At(kNtracksSector)))){
1216 AliWarning("No Histogram defined.");
1220 // TODO we should compare with
1221 // sector = Int_t(track->GetAlpha() / AliTRDgeometry::GetAlpha());
1223 AliTRDseedV1 *tracklet = NULL;
1225 for(Int_t itl = 0; itl < AliTRDgeometry::kNlayer; itl++){
1226 if(!(tracklet = fkTrack->GetTracklet(itl)) || !tracklet->IsOK()) continue;
1227 sector = static_cast<Int_t>(tracklet->GetDetector()/AliTRDgeometry::kNdets);
1235 //________________________________________________________
1236 void AliTRDcheckDET::SetRecoParam(AliTRDrecoParam *r)
1239 fReconstructor->SetRecoParam(r);
1242 //________________________________________________________
1243 void AliTRDcheckDET::GetDistanceToTracklet(Double_t *dist, AliTRDseedV1 * const tracklet, AliTRDcluster * const c)
1245 Float_t x = c->GetX();
1246 dist[0] = c->GetY() - tracklet->GetYat(x);
1247 dist[1] = c->GetZ() - tracklet->GetZat(x);
1250 //________________________________________________________
1251 void AliTRDcheckDET::GetEtaPhiAt(AliExternalTrackParam *track, Double_t x, Double_t &eta, Double_t &phi){
1253 // Get phi and eta at a given radial position
1255 AliExternalTrackParam workpar(*track);
1257 Double_t posLocal[3];
1258 Bool_t sucPos = workpar.GetXYZAt(x, fEventInfo->GetRunInfo()->GetMagneticField(), posLocal);
1259 Double_t sagPhi = sucPos ? TMath::ATan2(posLocal[1], posLocal[0]) : 0.;
1261 eta = workpar.Eta();
1265 //_______________________________________________________
1266 TH1* AliTRDcheckDET::MakePlotChi2()
1268 // Plot chi2/track normalized to number of degree of freedom
1269 // (tracklets) and compare with the theoretical distribution.
1271 // Alex Bercuci <A.Bercuci@gsi.de>
1275 TH2S *h2 = (TH2S*)fContainer->At(kChi2);
1276 TF1 f("fChi2", "[0]*pow(x, [1]-1)*exp(-0.5*x)", 0., 50.);
1277 f.SetParLimits(1,1, 1e100);
1278 TLegend *leg = new TLegend(.7,.7,.95,.95);
1279 leg->SetBorderSize(1); leg->SetHeader("Tracklets per Track");
1281 Bool_t kFIRST = kTRUE;
1282 for(Int_t il=1; il<=h2->GetNbinsX(); il++){
1283 h1 = h2->ProjectionY(Form("pyChi2%d", il), il, il);
1284 if(h1->Integral()<50) continue;
1285 h1->Scale(1./h1->Integral());
1286 h1->SetMarkerStyle(7);h1->SetMarkerColor(il);
1287 h1->SetLineColor(il);h1->SetLineStyle(2);
1288 f.SetParameter(1, .5*il);f.SetLineColor(il);
1289 h1->Fit(&f, "QW+", kFIRST ? "pc": "pcsame");
1290 leg->AddEntry(h1, Form("%d", il), "l");
1292 h1->GetXaxis()->SetRangeUser(0., 25.);
1302 //________________________________________________________
1303 TH1* AliTRDcheckDET::MakePlotNTracklets(){
1305 // Make nice bar plot of the number of tracklets in each method
1307 TH2F *tmp = (TH2F *)fContainer->FindObject("hNtlsBAR");
1308 TH1D *hBAR = tmp->ProjectionY();
1309 TH1F *hSTA = (TH1F *)fContainer->FindObject("hNtlsSTA");
1310 TH1F *hCON = (TH1F *)fContainer->FindObject("hNtls");
1311 TLegend *leg = new TLegend(0.13, 0.75, 0.39, 0.89);
1312 leg->SetBorderSize(1);
1313 leg->SetFillColor(0);
1315 Float_t scale = hCON->Integral();
1316 if(scale) hCON->Scale(100./scale);
1317 hCON->SetFillColor(kRed);hCON->SetLineColor(kRed);
1318 hCON->SetBarWidth(0.2);
1319 hCON->SetBarOffset(0.6);
1320 hCON->SetStats(kFALSE);
1321 hCON->GetYaxis()->SetRangeUser(0.,40.);
1322 hCON->GetYaxis()->SetTitleOffset(1.2);
1323 hCON->Draw("bar1"); leg->AddEntry(hCON, "Total", "f");
1324 hCON->SetMaximum(55.);
1326 if(scale) hBAR->Scale(100./scale);
1327 hBAR->SetFillColor(kGreen);hBAR->SetLineColor(kGreen);
1328 hBAR->SetBarWidth(0.2);
1329 hBAR->SetBarOffset(0.2);
1331 hBAR->SetStats(kFALSE);
1332 hBAR->GetYaxis()->SetRangeUser(0.,40.);
1333 hBAR->GetYaxis()->SetTitleOffset(1.2);
1334 hBAR->Draw("bar1same"); leg->AddEntry(hBAR, "Barrel", "f");
1336 if(scale) hSTA->Scale(100./scale);
1337 hSTA->SetFillColor(kBlue);hSTA->SetLineColor(kBlue);
1338 hSTA->SetBarWidth(0.2);
1339 hSTA->SetBarOffset(0.4);
1341 hSTA->SetStats(kFALSE);
1342 hSTA->GetYaxis()->SetRangeUser(0.,40.);
1343 hSTA->GetYaxis()->SetTitleOffset(1.2);
1344 hSTA->Draw("bar1same"); leg->AddEntry(hSTA, "Stand Alone", "f");
1350 //________________________________________________________
1351 void AliTRDcheckDET::MakePlotnTrackletsVsP(){
1353 // Plot abundance of tracks with number of tracklets as function of momentum
1359 Color_t color[AliTRDgeometry::kNlayer] = {kBlue, kOrange, kBlack, kGreen, kCyan, kRed};
1360 TH1 *h(NULL); TGraphErrors *g[AliTRDgeometry::kNlayer];
1361 for(Int_t itl(0); itl<AliTRDgeometry::kNlayer; itl++){
1362 g[itl] = new TGraphErrors();
1363 g[itl]->SetLineColor(color[itl]);
1364 g[itl]->SetMarkerColor(color[itl]);
1365 g[itl]->SetMarkerStyle(20 + itl);
1368 TH2 *hBar = (TH2F *)fContainer->FindObject("hNtlsBAR");
1369 TAxis *ax(hBar->GetXaxis());
1370 Int_t np(ax->GetNbins());
1371 for(Int_t ipBin(1); ipBin<np; ipBin++){
1372 h = hBar->ProjectionY("npBin", ipBin, ipBin);
1373 if(!Int_t(h->Integral())) continue;
1374 h->Scale(100./h->Integral());
1375 Float_t p(ax->GetBinCenter(ipBin));
1376 Float_t dp(ax->GetBinWidth(ipBin));
1377 Int_t ip(g[0]->GetN());
1378 for(Int_t itl(AliTRDgeometry::kNlayer); itl--;){
1379 g[itl]->SetPoint(ip, p, h->GetBinContent(itl+1));
1380 g[itl]->SetPointError(ip, dp/3.46, h->GetBinError(itl+1));
1384 TLegend *leg = new TLegend(0.76, 0.6, 1., 0.9);
1385 leg->SetBorderSize(0);
1386 leg->SetHeader("Tracklet/Track");
1387 leg->SetFillStyle(0);
1388 h = hBar->ProjectionX("npxBin"); h->Reset();
1390 h->GetYaxis()->SetRangeUser(1., 99.);
1391 h->GetYaxis()->SetMoreLogLabels();
1392 h->GetYaxis()->CenterTitle();
1393 h->GetYaxis()->SetTitleOffset(1.2);
1394 h->SetYTitle("Prob. [%]");
1395 h->GetXaxis()->SetRangeUser(0.4, 12.);
1396 h->GetXaxis()->SetMoreLogLabels();
1397 h->GetXaxis()->CenterTitle();
1399 for(Int_t itl(AliTRDgeometry::kNlayer); itl--;){
1401 leg->AddEntry(g[itl], Form("n = %d", itl+1),"pl");
1405 gPad->SetLogx();gPad->SetLogy();
1408 //________________________________________________________
1409 Bool_t AliTRDcheckDET::MakePlotPulseHeight(){
1411 // Create Plot of the Pluse Height Spectrum
1413 TCanvas *output = gPad->GetCanvas();
1417 TObjArray *arr = (TObjArray*)fContainer->FindObject("<PH>");
1418 h = (TH1F*)arr->At(0);
1419 h->SetMarkerStyle(24);
1420 h->SetMarkerColor(kBlack);
1421 h->SetLineColor(kBlack);
1422 h->GetYaxis()->SetTitleOffset(1.5);
1424 // Trending for the pulse height: plateau value, slope and timebin of the maximum
1425 TLinearFitter fit(1,"pol1");
1427 for(Int_t itime = 10; itime <= 20; itime++){
1428 time = static_cast<Double_t>(itime);
1429 fit.AddPoint(&time, h->GetBinContent(itime + 1), h->GetBinError(itime + 1));
1432 Double_t plateau = fit.GetParameter(0) + 12 * fit.GetParameter(1);
1433 Double_t slope = fit.GetParameter(1);
1434 PutTrendValue("PHplateau", plateau);
1435 PutTrendValue("PHslope", slope);
1436 PutTrendValue("PHamplificationPeak", static_cast<Double_t>(h->GetMaximumBin()-1));
1437 AliDebug(1, Form("plateau %f, slope %f, MaxTime %f", plateau, slope, static_cast<Double_t>(h->GetMaximumBin()-1)));
1438 // copy the second histogram in a new one with the same x-dimension as the phs with respect to time
1439 h1 = (TH1F *)arr->At(1);
1440 h2 = new TH1F("hphs1","Average PH", 31, -0.5, 30.5);
1441 for(Int_t ibin = h1->GetXaxis()->GetFirst(); ibin < h1->GetNbinsX(); ibin++)
1442 h2->SetBinContent(ibin, h1->GetBinContent(ibin));
1443 h2->SetMarkerStyle(22);
1444 h2->SetMarkerColor(kBlue);
1445 h2->SetLineColor(kBlue);
1448 // create axis according to the histogram dimensions of the original second histogram
1449 TGaxis *axis = new TGaxis(gPad->GetUxmin(),
1453 -0.08, 4.88, 510,"-L");
1454 axis->SetLineColor(kBlue);
1455 axis->SetLabelColor(kBlue);
1456 axis->SetTextColor(kBlue);
1457 axis->SetTitle("x_{0}-x_{c} [cm]");
1461 TH2 *ph2d = (TH2F *)arr->At(2);
1462 ph2d->GetYaxis()->SetTitleOffset(1.8);
1463 ph2d->SetStats(kFALSE);
1468 //________________________________________________________
1469 void AliTRDcheckDET::MakePlotMeanClustersLayer(){
1471 // Create Summary plot for the mean number of clusters per layer
1473 TCanvas *output = gPad->GetCanvas();
1474 output->Divide(3,2);
1475 TObjArray *histos = (TObjArray *)fContainer->At(kNclustersLayer);
1477 AliWarning("Histos for each layer not found");
1480 TProfile2D *hlayer = NULL;
1481 for(Int_t ily = 0; ily < AliTRDgeometry::kNlayer; ily++){
1482 hlayer = dynamic_cast<TProfile2D *>(histos->At(ily));
1483 output->cd(ily + 1);
1485 hlayer->Draw("colz");
1489 //________________________________________________________
1490 Bool_t AliTRDcheckDET::MakeBarPlot(TH1 *histo, Int_t color){
1492 // Draw nice bar plots
1494 if(!histo->GetEntries()) return kFALSE;
1495 histo->Scale(100./histo->Integral());
1496 histo->SetFillColor(color);
1497 histo->SetBarOffset(.2);
1498 histo->SetBarWidth(.6);
1499 histo->Draw("bar1");