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 ////////////////////////////////////////////////////////////////////////////
37 #include <TGraphErrors.h>
39 #include <TLinearFitter.h>
42 #include <TProfile2D.h>
43 #include <TObjArray.h>
45 #include <TObjString.h>
49 #include <TProfile2D.h>
54 #include "AliTRDcluster.h"
55 #include "AliESDHeader.h"
56 #include "AliESDRun.h"
57 #include "AliESDtrack.h"
58 #include "AliExternalTrackParam.h"
59 #include "AliTRDgeometry.h"
60 #include "AliTRDpadPlane.h"
61 #include "AliTRDSimParam.h"
62 #include "AliTRDseedV1.h"
63 #include "AliTRDtrackV1.h"
64 #include "AliTRDtrackerV1.h"
65 #include "AliTRDReconstructor.h"
66 #include "AliTrackReference.h"
67 #include "AliTrackPointArray.h"
68 #include "AliTracker.h"
69 #include "TTreeStream.h"
71 #include "info/AliTRDtrackInfo.h"
72 #include "info/AliTRDeventInfo.h"
73 #include "AliTRDinfoGen.h"
74 #include "AliTRDcheckDET.h"
79 ClassImp(AliTRDcheckDET)
81 //_______________________________________________________
82 AliTRDcheckDET::AliTRDcheckDET():
89 // Default constructor
91 SetNameTitle("TRDcheckDET", "Basic TRD data checker");
94 //_______________________________________________________
95 AliTRDcheckDET::AliTRDcheckDET(char* name):
96 AliTRDrecoTask(name, "Basic TRD data checker")
102 // Default constructor
108 //_______________________________________________________
109 AliTRDcheckDET::~AliTRDcheckDET(){
113 if(fTriggerNames) delete fTriggerNames;
116 //_______________________________________________________
117 void AliTRDcheckDET::UserCreateOutputObjects(){
119 // Create Output Objects
121 AliTRDrecoTask::UserCreateOutputObjects();
122 if(!fTriggerNames) fTriggerNames = new TMap();
125 //_______________________________________________________
126 void AliTRDcheckDET::UserExec(Option_t *opt){
128 // Execution function
129 // Filling TRD quality histos
131 AliDebug(2, Form("EventInfo[%p] Header[%p]", (void*)fEvent, (void*)(fEvent?fEvent->GetEventHeader():NULL)));
132 if(fEvent) fCentralityClass = fEvent->GetCentrality();
133 else fCentralityClass = -1; // Assume pp
135 AliTRDrecoTask::UserExec(opt);
137 TH1F *histo(NULL); AliTRDtrackInfo *fTrackInfo(NULL); Int_t nTracks(0); // Count the number of tracks per event
138 for(Int_t iti = 0; iti < fTracks->GetEntriesFast(); iti++){
139 if(!fTracks->UncheckedAt(iti)) continue;
140 if(!(fTrackInfo = dynamic_cast<AliTRDtrackInfo *>(fTracks->UncheckedAt(iti)))) continue;
141 if(!fTrackInfo->GetTrack()) continue;
145 if((histo = dynamic_cast<TH1F *>(fContainer->UncheckedAt(kNtracksEvent)))) histo->Fill(nTracks);
147 if(!fEvent->GetEventHeader()) return; // For trigger statistics event header is essential
148 Int_t triggermask = fEvent->GetEventHeader()->GetTriggerMask();
149 TString triggername = fEvent->GetRunInfo()->GetFiredTriggerClasses(triggermask);
150 AliDebug(6, Form("Trigger cluster: %d, Trigger class: %s\n", triggermask, triggername.Data()));
151 if((histo = dynamic_cast<TH1F *>(fContainer->UncheckedAt(kNeventsTrigger)))) histo->Fill(triggermask);
154 if((histo = dynamic_cast<TH1F *>(fContainer->UncheckedAt(kNeventsTriggerTracks)))) histo->Fill(triggermask);
156 if(triggermask <= 20 && !fTriggerNames->FindObject(Form("%d", triggermask))){
157 fTriggerNames->Add(new TObjString(Form("%d", triggermask)), new TObjString(triggername));
158 // also set the label for both histograms
159 if((histo = dynamic_cast<TH1F *>(fContainer->UncheckedAt(kNeventsTriggerTracks))))
160 histo->GetXaxis()->SetBinLabel(histo->FindBin(triggermask), triggername);
161 if((histo = dynamic_cast<TH1F *>(fContainer->UncheckedAt(kNeventsTrigger))))
162 histo->GetXaxis()->SetBinLabel(histo->FindBin(triggermask), triggername);
167 //_______________________________________________________
168 Bool_t AliTRDcheckDET::PostProcess(){
170 // Do Postprocessing (for the moment set the number of Reference histograms)
173 TH1F *h(NULL), *h1(NULL);
175 // Calculate of the trigger clusters purity
176 if((h = dynamic_cast<TH1F *>(fContainer->FindObject("hEventsTrigger"))) &&
177 (h1 = dynamic_cast<TH1F *>(fContainer->FindObject("hEventsTriggerTracks")))) {
179 Float_t purities[20], val = 0; memset(purities, 0, 20*sizeof(Float_t));
180 TString triggernames[20];
181 Int_t nTriggerClasses = 0;
182 for(Int_t ibin = 1; ibin <= h->GetNbinsX(); ibin++){
183 if((val = h1->GetBinContent(ibin))){
184 purities[nTriggerClasses] = val;
185 triggernames[nTriggerClasses] = h1->GetXaxis()->GetBinLabel(ibin);
190 if((h = dynamic_cast<TH1F *>(fContainer->UncheckedAt(kTriggerPurity)))){
191 TAxis *ax = h->GetXaxis();
192 for(Int_t itrg = 0; itrg < nTriggerClasses; itrg++){
193 h->Fill(itrg, purities[itrg]);
194 ax->SetBinLabel(itrg+1, triggernames[itrg].Data());
196 ax->SetRangeUser(-0.5, nTriggerClasses+.5);
197 h->GetYaxis()->SetRangeUser(0,1);
202 if((h=dynamic_cast<TH1F*>(fContainer->At(kTrackStatus)))){
203 Float_t ok = h->GetBinContent(1);
204 Int_t nerr = h->GetNbinsX();
205 for(Int_t ierr=nerr; ierr--;){
206 h->SetBinContent(ierr+1, ok>0.?1.e2*h->GetBinContent(ierr+1)/ok:0.);
208 h->SetBinContent(1, 0.);
212 TObjArray *arr(NULL);
213 if(( arr = dynamic_cast<TObjArray*>(fContainer->UncheckedAt(kTrackletStatus)))){
214 for(Int_t ily = AliTRDgeometry::kNlayer; ily--;){
215 if(!(h=dynamic_cast<TH1F*>(arr->At(ily)))) continue;
216 Float_t okB = h->GetBinContent(1);
217 Int_t nerrB = h->GetNbinsX();
218 for(Int_t ierr=nerrB; ierr--;){
219 h->SetBinContent(ierr+1, okB>0.?1.e2*h->GetBinContent(ierr+1)/okB:0.);
221 h->SetBinContent(1, 0.);
229 //_______________________________________________________
230 void AliTRDcheckDET::MakeSummary(){
232 // Create summary plots for TRD check DET
233 // This function creates 2 summary plots:
234 // - General Quantities
236 // The function will reuse GetRefFigure
239 TCanvas *cOut = new TCanvas(Form("summary%s1", GetName()), Form("Summary 1 for task %s", GetName()), 1024, 768);
242 // Create figures using GetRefFigure
243 cOut->cd(1); GetRefFigure(kFigNtracksEvent);
244 cOut->cd(2); GetRefFigure(kFigNtracksSector);
245 cOut->cd(3); GetRefFigure(kFigNclustersTrack);
246 cOut->cd(4); GetRefFigure(kFigNclustersTracklet);
247 cOut->cd(5); GetRefFigure(kFigNtrackletsTrack);
248 cOut->cd(6); GetRefFigure(kFigNTrackletsP);
249 cOut->cd(7); GetRefFigure(kFigChargeCluster);
250 cOut->cd(8); GetRefFigure(kFigChargeTracklet);
251 cOut->SaveAs("TRD_TrackQuality.gif");
255 cOut = new TCanvas(Form("summary%s2", GetName()), Form("Summary 2 for task %s", GetName()), 1024, 512);
256 cOut->cd(); GetRefFigure(kFigPH);
257 cOut->SaveAs("TRD_PH.gif");
260 // Third Plot: Mean Number of clusters as function of eta, phi and layer
261 cOut = new TCanvas(Form("summary%s3", GetName()), Form("Summary 3 for task %s", GetName()), 1024, 768);
262 cOut->cd(); MakePlotMeanClustersLayer();
263 cOut->SaveAs("TRD_MeanNclusters.gif");
268 //_______________________________________________________
269 Bool_t AliTRDcheckDET::GetRefFigure(Int_t ifig){
271 // Setting Reference Figures
275 TH1 *h = NULL; TObjArray *arr=NULL;
279 case kFigNclustersTrack:
280 (h=(TH1F*)fContainer->FindObject("hNcls"))->Draw("pl");
281 PutTrendValue("NClustersTrack", h->GetMean());
282 PutTrendValue("NClustersTrackRMS", h->GetRMS());
284 case kFigNclustersTracklet:
285 (h =(TH1F*)fContainer->FindObject("hNclTls"))->Draw("pc");
286 PutTrendValue("NClustersTracklet", h->GetMean());
287 PutTrendValue("NClustersTrackletRMS", h->GetRMS());
289 case kFigNtrackletsTrack:
290 h=MakePlotNTracklets();
292 PutTrendValue("NTrackletsTrack", h->GetMean());
293 PutTrendValue("NTrackletsTrackRMS", h->GetRMS());
296 case kFigNTrackletsP:
297 MakePlotnTrackletsVsP();
299 case kFigNtrackletsCross:
300 h = (TH1F*)fContainer->FindObject("hNtlsCross");
301 if(!MakeBarPlot(h, kRed)) break;
302 PutTrendValue("NTrackletsCross", h->GetMean());
303 PutTrendValue("NTrackletsCrossRMS", h->GetRMS());
305 case kFigNtrackletsFindable:
306 h = (TH1F*)fContainer->FindObject("hNtlsFindable");
307 if(!MakeBarPlot(h, kGreen)) break;
308 PutTrendValue("NTrackletsFindable", h->GetMean());
309 PutTrendValue("NTrackletsFindableRMS", h->GetRMS());
311 case kFigNtracksEvent:
312 (h = (TH1F*)fContainer->FindObject("hNtrks"))->Draw("pl");
313 PutTrendValue("NTracksEvent", h->GetMean());
314 PutTrendValue("NTracksEventRMS", h->GetRMS());
316 case kFigNtracksSector:
317 h = (TH1F*)fContainer->FindObject("hNtrksSector");
318 if(!MakeBarPlot(h, kGreen)) break;
319 PutTrendValue("NTracksSector", h->Integral()/h->GetNbinsX());
321 case kFigTrackStatus:
322 if(!(h=(TH1F *)fContainer->FindObject("hTrackStatus"))) break;
323 h->GetXaxis()->SetRangeUser(0.5, -1);
324 h->GetYaxis()->CenterTitle();
326 PutTrendValue("TrackStatus", h->Integral());
329 case kFigTrackletStatus:
330 if(!(arr = dynamic_cast<TObjArray*>(fContainer->At(kTrackletStatus)))) break;
331 leg = new TLegend(.68, .7, .97, .97);
332 leg->SetBorderSize(0);leg->SetFillStyle(0);
333 leg->SetHeader("TRD layer");
334 for(Int_t ily=AliTRDgeometry::kNlayer; ily--;){
335 if(!(h=dynamic_cast<TH1F*>(arr->At(ily)))) continue;
338 h->GetXaxis()->SetRangeUser(0.5, -1);
339 h->GetYaxis()->CenterTitle();
341 } else h->Draw("samepl");
342 leg->AddEntry(h, Form("ly = %d", ily), "l");
343 PutTrendValue(Form("TrackletStatus%d", ily), h->Integral());
352 gPad->SetMargin(0.125, 0.015, 0.1, 0.1);
353 MakePlotPulseHeight();
356 case kFigChargeCluster:
357 (h = (TH1F*)fContainer->FindObject("hQcl"))->Draw("c");
359 PutTrendValue("ChargeCluster", h->GetMaximumBin());
360 PutTrendValue("ChargeClusterRMS", h->GetRMS());
362 case kFigChargeTracklet:
363 (h=(TH1F*)fContainer->FindObject("hQtrklt"))->Draw("c");
364 PutTrendValue("ChargeTracklet", h->GetMaximumBin());
365 PutTrendValue("ChargeTrackletRMS", h->GetRMS());
367 case kFigNeventsTrigger:
368 ((TH1F*)fContainer->FindObject("hEventsTrigger"))->Draw("");
370 case kFigNeventsTriggerTracks:
371 ((TH1F*)fContainer->FindObject("hEventsTriggerTracks"))->Draw("");
373 case kFigTriggerPurity:
374 if(!MakeBarPlot((TH1F*)fContainer->FindObject("hTriggerPurity"), kGreen)) break;
379 AliInfo(Form("Reference plot [%d] missing result", ifig));
383 //_______________________________________________________
384 TObjArray *AliTRDcheckDET::Histos(){
386 // Create QA histograms
389 if(fContainer) return fContainer;
391 fContainer = new TObjArray(20);
392 //fContainer->SetOwner(kTRUE);
394 // Register Histograms
396 TH2 * h2 = NULL; // Pointer for two dimensional histograms
397 TH3 * h3 = NULL; // Pointer for tree dimensional histograms
399 if(!(h2 = (TH2F *)gROOT->FindObject("hNcls"))){
400 h2 = new TH2F("hNcls", "N_{clusters}/track;N_{clusters};Centrality;Entries", 181, -0.5, 180.5, AliTRDeventInfo::kCentralityClasses + 1, -1.5, AliTRDeventInfo::kCentralityClasses - 0.5);
402 fContainer->AddAt(h2, kNclustersTrack);
404 TObjArray *arr = new TObjArray(AliTRDgeometry::kNlayer);
405 arr->SetOwner(kTRUE); arr->SetName("clusters");
406 fContainer->AddAt(arr, kNclustersLayer);
407 for(Int_t ily=AliTRDgeometry::kNlayer; ily--;){
408 if(!(h = (TProfile2D *)gROOT->FindObject(Form("hNcl%d", ily)))){
409 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());
410 h->GetXaxis()->SetTitle("#eta");
411 h->GetYaxis()->SetTitle("#phi");
416 if(!(h2 = (TH2F *)gROOT->FindObject("hNclTls"))){
417 h2 = new TH2F("hNclTls","N_{clusters}/tracklet;N_{clusters};Entries", 51, -0.5, 50.5, AliTRDeventInfo::kCentralityClasses + 1, -1.5, AliTRDeventInfo::kCentralityClasses - 0.5);
419 fContainer->AddAt(h2, kNclustersTracklet);
421 if(!(h2 = (TH2F *)gROOT->FindObject("hNtls"))){
422 h2 = new TH2F("hNtls", "N_{tracklets}/track;N^{tracklet};Centrality;freq.[%]", AliTRDgeometry::kNlayer, 0.5, 6.5, AliTRDeventInfo::kCentralityClasses + 1, -1.5, AliTRDeventInfo::kCentralityClasses - 0.5);
424 fContainer->AddAt(h2, kNtrackletsTrack);
426 if(!(h = (TH2F *)gROOT->FindObject("htlsSTA"))){
427 h = new TH2F("hNtlsSTA", "#splitline{N_{tracklets}/track}{Stand Alone};N^{tracklet};Centrality;freq.[%]", AliTRDgeometry::kNlayer, 0.5, 6.5, AliTRDeventInfo::kCentralityClasses + 1, -1.5, AliTRDeventInfo::kCentralityClasses - 0.5);
429 fContainer->AddAt(h, kNtrackletsSTA);
431 // Binning for momentum dependent tracklet Plots
434 Float_t binsP[kNp+1], binsTrklt[AliTRDgeometry::kNlayer+1], binsCent[AliTRDeventInfo::kCentralityClasses+2];
435 for(Int_t i=0;i<kNp+1; i++,p+=(TMath::Exp(i*i*.001)-1.)) binsP[i]=p;
436 for(Int_t il = 0; il <= AliTRDgeometry::kNlayer; il++) binsTrklt[il] = 0.5 + il;
437 for(Int_t icent = -1; icent < AliTRDeventInfo::kCentralityClasses + 1; icent++) binsCent[icent+1] = icent - 0.5;
438 if(!(h3 = (TH3F *)gROOT->FindObject("htlsBAR"))){
439 // Make tracklets for barrel tracks momentum dependent (if we do not exceed min and max values)
440 h3 = new TH3F("hNtlsBAR",
441 "N_{tracklets}/track;p [GeV/c];N^{tracklet};freq. [%]",
442 kNp, binsP, AliTRDgeometry::kNlayer, binsTrklt, AliTRDeventInfo::kCentralityClasses + 1, binsCent);
444 fContainer->AddAt(h3, kNtrackletsBAR);
447 if(!(h2 = (TH2F *)gROOT->FindObject("hNtlsCross"))){
448 h2 = new TH2F("hNtlsCross", "N_{tracklets}^{cross}/track;n_{row cross};Centrality;freq.[%]", 7, -0.5, 6.5, AliTRDeventInfo::kCentralityClasses + 1, -1.5, AliTRDeventInfo::kCentralityClasses - 0.5);
450 fContainer->AddAt(h2, kNtrackletsCross);
452 if(!(h2 = (TH2F *)gROOT->FindObject("hNtlsFindable"))){
453 h2 = new TH2F("hNtlsFindable", "Found/Findable Tracklets;r[a.u];Centrality;Entries" , 101, -0.005, 1.005, AliTRDeventInfo::kCentralityClasses + 1, -1.5, AliTRDeventInfo::kCentralityClasses - 0.5);
455 fContainer->AddAt(h2, kNtrackletsFindable);
457 if(!(h = (TH1F *)gROOT->FindObject("hNtrks"))){
458 h = new TH1F("hNtrks", "N_{tracks} / event", 100, 0, 100);
459 h->GetXaxis()->SetTitle("N_{tracks}");
460 h->GetYaxis()->SetTitle("Entries");
462 fContainer->AddAt(h, kNtracksEvent);
464 if(!(h = (TH1F *)gROOT->FindObject("hNtrksSector"))){
465 h = new TH1F("hNtrksSector", "N_{tracks} / sector", AliTRDgeometry::kNsector, -0.5, 17.5);
466 h->GetXaxis()->SetTitle("sector");
467 h->GetYaxis()->SetTitle("freq. [%]");
469 fContainer->AddAt(h, kNtracksSector);
471 if(!(h = (TH1F*)gROOT->FindObject("hTrackStatus"))){
472 const Int_t nerr = 7;
473 h = new TH1F("hTrackStatus", "Track Status", nerr, -0.5, nerr-0.5);
474 const Char_t *label[nerr] = {"OK", "PROL", "PROP", "AJST", "SNP", "TINI", "UPDT"};
476 for(Int_t ierr=nerr; ierr--;) ax->SetBinLabel(ierr+1, label[ierr]);
477 h->SetYTitle("Relative Error to Good [%]");
479 fContainer->AddAt(h, kTrackStatus);
481 arr = new TObjArray(AliTRDgeometry::kNlayer);
482 arr->SetOwner(kTRUE); arr->SetName("TrackletStatus");
483 fContainer->AddAt(arr, kTrackletStatus);
484 for(Int_t ily=AliTRDgeometry::kNlayer; ily--;){
485 if(!(h = (TH1F *)gROOT->FindObject(Form("hTrackletStatus%d", ily)))){
486 const Int_t nerr = 8;
487 h = new TH1F(Form("hTrackletStatus%d", ily), "Tracklet status", nerr, -0.5, nerr-0.5);
488 h->SetLineColor(ily+1);
489 const Char_t *label[nerr] = {"OK", "Geom", "Bound", "NoCl", "NoAttach", "NoClTr", "NoFit", "Chi2"};
491 for(Int_t ierr=nerr; ierr--;) ax->SetBinLabel(ierr+1, label[ierr]);
492 h->SetYTitle("Relative Error to Good [%]");
498 arr = new TObjArray(3);
499 arr->SetOwner(kTRUE); arr->SetName("<PH>");
500 fContainer->AddAt(arr, kPH);
501 if(!(h3 = (TH3F *)gROOT->FindObject("hPHx"))){
502 h3 = new TH3F("hPHx", "<PH>(x);x[mm];Centrality;Charge[a.u.]", 31, -0.08, 4.88, 100, 0, 1024, AliTRDeventInfo::kCentralityClasses + 1, -1.5, AliTRDeventInfo::kCentralityClasses - 0.5);
505 if(!(h3 = (TH3F *)gROOT->FindObject("hPHt"))){
506 h3 = new TH3F("hPHt", "<PH>(t);time[100ns];Centrality;Charge[a.u.]", 31, -0.5, 30.5, 100, 0, 1024, AliTRDeventInfo::kCentralityClasses + 1, -1.5, AliTRDeventInfo::kCentralityClasses - 0.5);
511 if(!(h3 = (TH3S*)gROOT->FindObject("hChi2"))){
512 h3 = new TH3S("hChi2", "#chi^{2}/track;ndf;#chi^{2}/ndf;Centrality", AliTRDgeometry::kNlayer, .5, AliTRDgeometry::kNlayer+.5, 100, 0, 50, AliTRDeventInfo::kCentralityClasses + 1, -1.5, AliTRDeventInfo::kCentralityClasses - 0.5);
514 fContainer->AddAt(h3, kChi2);
516 if(!(h2 = (TH2F *)gROOT->FindObject("hQcl"))){
517 h2 = new TH2F("hQcl", "Q_{cluster};Charge[a.u.];Centrality;Entries", 200, 0, 1200, AliTRDeventInfo::kCentralityClasses + 1, -1.5, AliTRDeventInfo::kCentralityClasses - 0.5);
519 fContainer->AddAt(h2, kChargeCluster);
521 if(!(h2 = (TH2F *)gROOT->FindObject("hQtrklt"))){
522 h2 = new TH2F("hQtrklt", "Q_{tracklet};Charge[a.u.];Centrality;Entries", 6000, 0, 6000, AliTRDeventInfo::kCentralityClasses + 1, -1.5, AliTRDeventInfo::kCentralityClasses - 0.5);
524 fContainer->AddAt(h2, kChargeTracklet);
527 if(!(h = (TH1F *)gROOT->FindObject("hEventsTrigger")))
528 h = new TH1F("hEventsTrigger", "Trigger Class", 100, 0, 100);
531 fContainer->AddAt(h, kNeventsTrigger);
533 if(!(h = (TH1F *)gROOT->FindObject("hEventsTriggerTracks")))
534 h = new TH1F("hEventsTriggerTracks", "Trigger Class (Tracks)", 100, 0, 100);
536 fContainer->AddAt(h, kNeventsTriggerTracks);
538 if(!(h = (TH1F *)gROOT->FindObject("hTriggerPurity"))){
539 h = new TH1F("hTriggerPurity", "Trigger Purity", 10, -0.5, 9.5);
540 h->GetXaxis()->SetTitle("Trigger Cluster");
541 h->GetYaxis()->SetTitle("freq.");
543 fContainer->AddAt(h, kTriggerPurity);
552 //_______________________________________________________
553 TH1 *AliTRDcheckDET::PlotTrackStatus(const AliTRDtrackV1 *track)
556 // Plot the track propagation status. The following errors are defined (see AliTRDtrackV1::ETRDtrackError)
557 // PROL - track prolongation failure
558 // PROP - track propagation failure
559 // AJST - crossing sectors failure
560 // SNP - too large bending
561 // TINI - tracklet initialization failure
562 // UPDT - track position/covariance update failure
564 // Performance plot looks as below:
566 //<img src="TRD/trackStatus.gif">
569 if(track) fkTrack = track;
571 AliDebug(4, "No Track defined.");
575 if(!(h = dynamic_cast<TH1F *>(fContainer->At(kTrackStatus)))){
576 AliWarning("No Histogram defined.");
579 h->Fill(fkTrack->GetStatusTRD());
583 //_______________________________________________________
584 TH1 *AliTRDcheckDET::PlotTrackletStatus(const AliTRDtrackV1 *track)
587 // Plot the tracklet propagation status. The following errors are defined for tracklet (see AliTRDtrackV1::ETRDlayerError)
589 // Bound - tracklet too close to chamber walls
590 // NoCl - no clusters in the track roads
591 // NoAttach - fail to attach clusters
592 // NoClTr - fail to use clusters for fit
593 // NoFit - tracklet fit failled
594 // Chi2 - chi2 tracklet-track over threshold
596 // Performance plot looks as below:
598 //<img src="TRD/trackletStatus.gif">
601 if(track) fkTrack = track;
603 AliDebug(4, "No Track defined.");
606 TObjArray *arr =NULL;
607 if(!(arr = dynamic_cast<TObjArray*>(fContainer->At(kTrackletStatus)))){
608 AliWarning("Histograms not defined.");
613 for(Int_t ily=AliTRDgeometry::kNlayer; ily--;){
614 if(!(h = dynamic_cast<TH1F*>(arr->At(ily)))){
615 AliWarning(Form("Missing histo for layer %d.", ily));
618 h->Fill(fkTrack->GetStatusTRD(ily));
623 //_______________________________________________________
624 TH1 *AliTRDcheckDET::PlotNClustersTracklet(const AliTRDtrackV1 *track){
626 // Plot the mean number of clusters per tracklet
628 if(track) fkTrack = track;
630 AliDebug(4, "No Track defined.");
633 AliExternalTrackParam *par = fkTrack->GetTrackIn() ? fkTrack->GetTrackIn() : fkTrack->GetTrackOut();
635 TProfile2D *hlayer = NULL;
636 Double_t eta = 0., phi = 0.;
637 if(!(h = dynamic_cast<TH2F *>(fContainer->At(kNclustersTracklet)))){
638 AliWarning("No Histogram defined.");
641 AliTRDseedV1 *tracklet = NULL;
642 TObjArray *histosLayer = dynamic_cast<TObjArray *>(fContainer->At(kNclustersLayer));
644 AliWarning("No Histograms for single layer defined");
646 for(Int_t itl = 0; itl < AliTRDgeometry::kNlayer; itl++){
647 if(!(tracklet = fkTrack->GetTracklet(itl)) || !tracklet->IsOK()) continue;
648 h->Fill(tracklet->GetN2(), fCentralityClass);
649 if(histosLayer && par){
650 if((hlayer = dynamic_cast<TProfile2D *>(histosLayer->At(itl)))){
651 GetEtaPhiAt(par, tracklet->GetX0(), eta, phi);
652 hlayer->Fill(eta, phi, tracklet->GetN2());
659 //_______________________________________________________
660 TH1 *AliTRDcheckDET::PlotNClustersTrack(const AliTRDtrackV1 *track){
662 // Plot the number of clusters in one track
664 if(track) fkTrack = track;
666 AliDebug(4, "No Track defined.");
670 if(!(h = dynamic_cast<TH2F *>(fContainer->At(kNclustersTrack)))){
671 AliWarning("No Histogram defined.");
676 AliTRDseedV1 *tracklet = NULL;
677 AliExternalTrackParam *par = fkTrack->GetTrackOut() ? fkTrack->GetTrackOut() : fkTrack->GetTrackIn();
678 if(!par) return NULL;
679 Double_t momentumRec = par->P();
680 for(Int_t itl = 0; itl < AliTRDgeometry::kNlayer; itl++){
681 if(!(tracklet = fkTrack->GetTracklet(itl)) || !tracklet->IsOK()) continue;
682 Int_t n(tracklet->GetN());
684 if(DebugLevel() > 2){
685 Int_t crossing = Int_t(tracklet->IsRowCross());
686 Int_t detector = tracklet->GetDetector();
687 Float_t theta = TMath::ATan(tracklet->GetZref(1));
688 Float_t phi = TMath::ATan(tracklet->GetYref(1));
689 Float_t momentumMC = 0.;
691 Int_t kinkIndex = fkESD ? fkESD->GetKinkIndex() : 0;
692 UShort_t nclsTPC = fkESD ? fkESD->GetTPCncls() : 0;
694 if(fkMC->GetTrackRef()) momentumMC = fkMC->GetTrackRef()->P();
695 pdg = fkMC->GetPDG();
697 (*DebugStream()) << "NClustersTrack"
698 << "Detector=" << detector
699 << "Centrality="<< fCentralityClass
700 << "crossing=" << crossing
701 << "momentumMC="<< momentumMC
702 << "momentumRec="<< momentumRec
706 << "kinkIndex=" << kinkIndex
707 << "TPCncls=" << nclsTPC
712 h->Fill(nclusters, fCentralityClass);
717 //_______________________________________________________
718 TH1 *AliTRDcheckDET::PlotNTrackletsTrack(const AliTRDtrackV1 *track){
720 // Plot the number of tracklets
722 if(track) fkTrack = track;
724 AliDebug(4, "No Track defined.");
727 TH2 *h = NULL, *hSta = NULL; TH3 *hBarrel = NULL;
728 if(!(h = dynamic_cast<TH2F *>(fContainer->At(kNtrackletsTrack)))){
729 AliWarning("No Histogram defined.");
732 Int_t nTracklets = fkTrack->GetNumberOfTracklets();
733 h->Fill(nTracklets, fCentralityClass);
735 Int_t status = fkESD->GetStatus();
737 /* 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);*/
739 Int_t method = -1; // to distinguish between stand alone and full barrel tracks in the debugging
740 if((status & AliESDtrack::kTRDin) != 0){
742 // Full Barrel Track: Save momentum dependence
743 if(!(hBarrel = dynamic_cast<TH3F *>(fContainer->At(kNtrackletsBAR)))){
744 AliWarning("Method: Barrel. Histogram not processed!");
747 AliExternalTrackParam *par(fkTrack->GetTrackIn());
749 AliError("Input track params missing");
752 p = par->P(); // p needed later in the debug streaming
753 hBarrel->Fill(p, nTracklets, fCentralityClass);
755 // Stand alone Track: momentum dependence not usefull
757 if(!(hSta = dynamic_cast<TH2F *>(fContainer->At(kNtrackletsSTA)))) {
758 AliWarning("Method: StandAlone. Histogram not processed!");
761 hSta->Fill(nTracklets, fCentralityClass);
764 if(DebugLevel() > 2){
765 AliTRDseedV1 *tracklet = NULL;
766 AliTRDgeometry *geo(AliTRDinfoGen::Geometry());
767 Int_t sector = -1, stack = -1, detector;
768 for(Int_t itl = 0; itl < AliTRDgeometry::kNlayer; itl++){
769 if(!(tracklet = fkTrack->GetTracklet(itl)) || !(tracklet->IsOK())) continue;
770 detector = tracklet->GetDetector();
771 sector = geo->GetSector(detector);
772 stack = geo->GetStack(detector);
775 (*DebugStream()) << "NTrackletsTrack"
776 << "Sector=" << sector
778 << "Centrality=" << fCentralityClass
779 << "NTracklets=" << nTracklets
780 << "Method=" << method
784 if(DebugLevel() > 3){
785 AliTRDseedV1 *tracklet = NULL;
786 for(Int_t il = 0; il < AliTRDgeometry::kNlayer; il++){
787 if((tracklet = fkTrack->GetTracklet(il)) && tracklet->IsOK()){
788 (*DebugStream()) << "NTrackletsLayer"
799 //_______________________________________________________
800 TH1 *AliTRDcheckDET::PlotNTrackletsRowCross(const AliTRDtrackV1 *track){
802 // Plot the number of tracklets
804 if(track) fkTrack = track;
806 AliDebug(4, "No Track defined.");
810 if(!(h = dynamic_cast<TH2F *>(fContainer->At(kNtrackletsCross)))){
811 AliWarning("No Histogram defined.");
816 AliTRDseedV1 *tracklet = NULL;
817 for(Int_t il = 0; il < AliTRDgeometry::kNlayer; il++){
818 if(!(tracklet = fkTrack->GetTracklet(il)) || !tracklet->IsOK()) continue;
820 if(tracklet->IsRowCross()) ncross++;
822 h->Fill(ncross, fCentralityClass);
826 //_______________________________________________________
827 TH1 *AliTRDcheckDET::PlotFindableTracklets(const AliTRDtrackV1 *track){
829 // Plots the ratio of number of tracklets vs.
830 // number of findable tracklets
832 // Findable tracklets are defined as track prolongation
833 // to layer i does not hit the dead area +- epsilon
835 // In order to check whether tracklet hist active area in Layer i,
836 // the track is refitted and the fitted position + an uncertainty
837 // range is compared to the chamber border (also with a different
840 // For the track fit two cases are distinguished:
841 // If the track is a stand alone track (defined by the status bit
842 // encoding, then the track is fitted with the tilted Rieman model
843 // Otherwise the track is fitted with the Kalman fitter in two steps:
844 // Since the track parameters are give at the outer point, we first
845 // fit in direction inwards. Afterwards we fit again in direction outwards
846 // to extrapolate the track to layers which are not reached by the track
847 // For the Kalman model, the radial track points have to be shifted by
848 // a distance epsilon in the direction that we want to fit
850 const Float_t epsilon = 0.01; // dead area tolerance
851 const Float_t epsilonR = 1; // shift in radial direction of the anode wire position (Kalman filter only)
852 const Float_t deltaY = 0.7; // Tolerance in the track position in y-direction
853 const Float_t deltaZ = 7.0; // Tolerance in the track position in z-direction (Padlength)
854 Double_t xAnode[AliTRDgeometry::kNlayer] = {300.2, 312.8, 325.4, 338.0, 350.6, 363.2}; // Take the default X0
856 if(track) fkTrack = track;
858 AliDebug(4, "No Track defined.");
862 if(!(h = dynamic_cast<TH2F *>(fContainer->At(kNtrackletsFindable)))){
863 AliWarning("No Histogram defined.");
866 Int_t nFound = 0, nFindable = 0;
868 Double_t ymin = 0., ymax = 0., zmin = 0., zmax = 0.;
869 Double_t y = 0., z = 0.;
870 AliTRDseedV1 *tracklet = NULL;
872 for(Int_t il = 0; il < AliTRDgeometry::kNlayer; il++){
873 if((tracklet = fkTrack->GetTracklet(il)) && tracklet->IsOK()){
874 tracklet->SetReconstructor(AliTRDinfoGen::Reconstructor());
878 // 2 Different cases:
879 // 1st stand alone: here we cannot propagate, but be can do a Tilted Rieman Fit
880 // 2nd barrel track: here we propagate the track to the layers
881 AliTrackPoint points[6];
883 memset(xyz, 0, sizeof(Float_t) * 3);
884 if(((fkESD->GetStatus() & AliESDtrack::kTRDout) > 0) && !((fkESD->GetStatus() & AliESDtrack::kTRDin) > 0)){
886 for(Int_t il = 0; il < AliTRDgeometry::kNlayer; il++){
888 points[il].SetXYZ(xyz);
890 AliTRDtrackerV1::FitRiemanTilt(const_cast<AliTRDtrackV1 *>(fkTrack), NULL, kTRUE, 6, points);
896 // -> Kalman outwards
897 AliTRDtrackV1 copyTrack(*fkTrack); // Do Kalman on a (non-constant) copy of the track
898 AliTrackPoint pointsInward[6], pointsOutward[6];
899 for(Int_t il = AliTRDgeometry::kNlayer; il--;){
900 // In order to avoid complications in the Kalman filter if the track points have the same radial
901 // position like the tracklets, we have to shift the radial postion of the anode wire by epsilon
902 // in the direction we want to go
903 // The track points have to be in reverse order for the Kalman Filter inwards
904 xyz[0] = xAnode[AliTRDgeometry::kNlayer - il - 1] - epsilonR;
905 pointsInward[il].SetXYZ(xyz);
906 xyz[0] = xAnode[il] + epsilonR;
907 pointsOutward[il].SetXYZ(xyz);
909 /*for(Int_t ipt = 0; ipt < AliTRDgeometry::kNlayer; ipt++)
910 printf("%d. X = %f\n", ipt, points[ipt].GetX());*/
912 AliTRDtrackerV1::FitKalman(©Track, NULL, kFALSE, 6, pointsInward);
913 memcpy(points, pointsInward, sizeof(AliTrackPoint) * 6); // Preliminary store the inward results in the Array points
915 AliTRDtrackerV1::FitKalman(©Track, NULL, kTRUE, 6, pointsInward);
916 memcpy(points, pointsOutward, sizeof(AliTrackPoint) * AliTRDgeometry::kNlayer);
918 AliTRDgeometry *geo(AliTRDinfoGen::Geometry());
919 for(Int_t il = 0; il < AliTRDgeometry::kNlayer; il++){
920 y = points[il].GetY();
921 z = points[il].GetZ();
922 if((stack = geo->GetStack(z, il)) < 0) continue; // Not findable
923 pp = geo->GetPadPlane(il, stack);
924 ymin = pp->GetCol0() + epsilon;
925 ymax = pp->GetColEnd() - epsilon;
926 zmin = pp->GetRowEnd() + epsilon;
927 zmax = pp->GetRow0() - epsilon;
928 // ignore y-crossing (material)
929 if((z + deltaZ > zmin && z - deltaZ < zmax) && (y + deltaY > ymin && y - deltaY < ymax)) nFindable++;
930 if(DebugLevel() > 3){
931 Double_t posTracklet[2] = {tracklet ? tracklet->GetYfit(0) : 0, tracklet ? tracklet->GetZfit(0) : 0};
932 Int_t hasTracklet = tracklet ? 1 : 0;
933 (*DebugStream()) << "FindableTracklets"
935 << "ytracklet=" << posTracklet[0]
937 << "ztracklet=" << posTracklet[1]
939 << "tracklet=" << hasTracklet
944 h->Fill((nFindable > 0 ? TMath::Min(nFound/static_cast<Double_t>(nFindable), 1.) : 1), fCentralityClass);
945 AliDebug(2, Form("Findable[Found]: %d[%d|%f]", nFindable, nFound, nFound/static_cast<Float_t>(nFindable > 0 ? nFindable : 1)));
950 //_______________________________________________________
951 TH1 *AliTRDcheckDET::PlotChi2(const AliTRDtrackV1 *track){
953 // Plot the chi2 of the track
955 if(track) fkTrack = track;
957 AliDebug(4, "No Track defined.");
961 if(!(h = dynamic_cast<TH3S*>(fContainer->At(kChi2)))) {
962 AliWarning("No Histogram defined.");
965 Int_t n = fkTrack->GetNumberOfTracklets();
968 h->Fill(n, fkTrack->GetChi2()/n, fCentralityClass);
973 //_______________________________________________________
974 TH1 *AliTRDcheckDET::PlotPHt(const AliTRDtrackV1 *track){
976 // Plot the average pulse height
978 if(track) fkTrack = track;
980 AliDebug(4, "No Track defined.");
984 if(!(h = dynamic_cast<TH3F *>(((TObjArray*)(fContainer->At(kPH)))->At(1)))){
985 AliWarning("No Histogram defined.");
988 AliTRDseedV1 *tracklet = NULL;
989 AliTRDcluster *c = NULL;
990 for(Int_t itl = 0; itl < AliTRDgeometry::kNlayer; itl++){
991 if(!(tracklet = fkTrack->GetTracklet(itl)) || !tracklet->IsOK())continue;
992 Int_t crossing = Int_t(tracklet->IsRowCross());
993 Int_t detector = tracklet->GetDetector();
994 tracklet->ResetClusterIter();
995 while((c = tracklet->NextCluster())){
996 if(!IsUsingClustersOutsideChamber() && !c->IsInChamber()) continue;
997 Int_t localtime = c->GetLocalTimeBin();
998 Double_t absoluteCharge = TMath::Abs(c->GetQ());
999 h->Fill(localtime, absoluteCharge, fCentralityClass);
1000 if(DebugLevel() > 3){
1001 Int_t inChamber = c->IsInChamber() ? 1 : 0;
1002 Double_t distance[2];
1003 GetDistanceToTracklet(distance, tracklet, c);
1004 Float_t theta = TMath::ATan(tracklet->GetZref(1));
1005 Float_t phi = TMath::ATan(tracklet->GetYref(1));
1006 AliExternalTrackParam *trdPar = fkTrack->GetTrackIn();
1007 Float_t momentumMC = 0, momentumRec = trdPar ? trdPar->P() : fkTrack->P(); // prefer Track Low
1009 Int_t kinkIndex = fkESD ? fkESD->GetKinkIndex() : 0;
1010 UShort_t tpcNCLS = fkESD ? fkESD->GetTPCncls() : 0;
1012 if(fkMC->GetTrackRef()) momentumMC = fkMC->GetTrackRef()->P();
1013 pdg = fkMC->GetPDG();
1015 (*DebugStream()) << "PHt"
1016 << "Detector=" << detector
1017 << "Centrality="<< fCentralityClass
1018 << "crossing=" << crossing
1019 << "inChamber=" << inChamber
1020 << "Timebin=" << localtime
1021 << "Charge=" << absoluteCharge
1022 << "momentumMC=" << momentumMC
1023 << "momentumRec=" << momentumRec
1025 << "theta=" << theta
1027 << "kinkIndex=" << kinkIndex
1028 << "TPCncls=" << tpcNCLS
1029 << "dy=" << distance[0]
1030 << "dz=" << distance[1]
1039 //_______________________________________________________
1040 TH1 *AliTRDcheckDET::PlotPHx(const AliTRDtrackV1 *track){
1042 // Plots the average pulse height vs the distance from the anode wire
1043 // (plus const anode wire offset)
1045 if(track) fkTrack = track;
1047 AliDebug(4, "No Track defined.");
1051 if(!(h = dynamic_cast<TH3F *>(((TObjArray*)(fContainer->At(kPH)))->At(0)))){
1052 AliWarning("No Histogram defined.");
1055 AliTRDseedV1 *tracklet(NULL);
1056 AliTRDcluster *c(NULL);
1057 Double_t xd(0.), dqdl(0.);
1058 TVectorD vq(AliTRDseedV1::kNtb), vxd(AliTRDseedV1::kNtb), vdqdl(AliTRDseedV1::kNtb);
1059 for(Int_t itl = 0; itl < AliTRDgeometry::kNlayer; itl++){
1060 if(!(tracklet = fkTrack->GetTracklet(itl)) || !(tracklet->IsOK())) continue;
1061 Int_t det(tracklet->GetDetector());
1062 Bool_t rc(tracklet->IsRowCross());
1063 for(Int_t ic(0); ic<AliTRDseedV1::kNtb; ic++){
1064 Bool_t kFIRST(kFALSE);
1065 if(!(c = tracklet->GetClusters(ic))){
1066 if(!(c = tracklet->GetClusters(AliTRDseedV1::kNtb+ic))) continue;
1067 } else kFIRST=kTRUE;
1068 if(!IsUsingClustersOutsideChamber() && !c->IsInChamber()) continue;
1069 xd = tracklet->GetX0() - c->GetX(); vxd[ic] = xd;
1070 dqdl=tracklet->GetdQdl(ic); vdqdl[ic] = dqdl;
1072 if(kFIRST && (c = tracklet->GetClusters(AliTRDseedV1::kNtb+ic))) vq[ic]+=c->GetQ();
1073 h->Fill(xd, dqdl, fCentralityClass);
1075 if(DebugLevel() > 3){
1076 (*DebugStream()) << "PHx"
1081 << "dqdl=" << &vdqdl
1088 //_______________________________________________________
1089 TH1 *AliTRDcheckDET::PlotChargeCluster(const AliTRDtrackV1 *track){
1091 // Plot the cluster charge
1093 if(track) fkTrack = track;
1095 AliDebug(4, "No Track defined.");
1099 if(!(h = dynamic_cast<TH2F *>(fContainer->At(kChargeCluster)))){
1100 AliWarning("No Histogram defined.");
1103 AliTRDseedV1 *tracklet = NULL;
1104 AliTRDcluster *c = NULL;
1105 for(Int_t itl = 0; itl < AliTRDgeometry::kNlayer; itl++){
1106 if(!(tracklet = fkTrack->GetTracklet(itl)) || !tracklet->IsOK())continue;
1107 for(Int_t ic(0); ic < AliTRDseedV1::kNtb; ic++){
1108 Bool_t kFIRST(kFALSE);
1109 if(!(c = tracklet->GetClusters(ic))) {
1110 if(!(c = tracklet->GetClusters(AliTRDseedV1::kNtb+ic))) continue;
1111 } else kFIRST = kTRUE;
1112 Float_t q(c->GetQ());
1113 if(kFIRST && (c = tracklet->GetClusters(AliTRDseedV1::kNtb+ic))) q+=c->GetQ();
1114 h->Fill(q, fCentralityClass);
1120 //_______________________________________________________
1121 TH1 *AliTRDcheckDET::PlotChargeTracklet(const AliTRDtrackV1 *track){
1123 // Plot the charge deposit per chamber
1125 if(track) fkTrack = track;
1127 AliDebug(4, "No Track defined.");
1131 if(!(h = dynamic_cast<TH2F *>(fContainer->At(kChargeTracklet)))){
1132 AliWarning("No Histogram defined.");
1135 AliTRDseedV1 *tracklet = NULL;
1136 AliTRDcluster *c = NULL;
1138 Int_t nTracklets =fkTrack->GetNumberOfTracklets();
1139 for(Int_t itl(0); itl < AliTRDgeometry::kNlayer; itl++){
1140 if(!(tracklet = fkTrack->GetTracklet(itl)) || !tracklet->IsOK()) continue;
1142 for(Int_t ic = AliTRDseedV1::kNclusters; ic--;){
1143 if(!(c = tracklet->GetClusters(ic))) continue;
1144 qTot += TMath::Abs(c->GetQ());
1146 h->Fill(qTot, fCentralityClass);
1147 if(DebugLevel() > 3){
1148 Int_t crossing = (Int_t)tracklet->IsRowCross();
1149 Int_t detector = tracklet->GetDetector();
1150 Float_t theta = TMath::ATan(tracklet->GetZfit(1));
1151 Float_t phi = TMath::ATan(tracklet->GetYfit(1));
1152 Float_t momentum = 0.;
1154 Int_t kinkIndex = fkESD ? fkESD->GetKinkIndex() : 0;
1155 UShort_t nclsTPC = fkESD ? fkESD->GetTPCncls() : 0;
1157 if(fkMC->GetTrackRef()) momentum = fkMC->GetTrackRef()->P();
1158 pdg = fkMC->GetPDG();
1160 (*DebugStream()) << "ChargeTracklet"
1161 << "Detector=" << detector
1162 << "Centrality="<< fCentralityClass
1163 << "crossing=" << crossing
1164 << "momentum=" << momentum
1165 << "nTracklets="<< nTracklets
1167 << "theta=" << theta
1169 << "kinkIndex=" << kinkIndex
1170 << "TPCncls=" << nclsTPC
1178 //_______________________________________________________
1179 TH1 *AliTRDcheckDET::PlotNTracksSector(const AliTRDtrackV1 *track){
1181 // Plot the number of tracks per Sector
1183 if(track) fkTrack = track;
1185 AliDebug(4, "No Track defined.");
1189 if(!(h = dynamic_cast<TH1F *>(fContainer->At(kNtracksSector)))){
1190 AliWarning("No Histogram defined.");
1194 // TODO we should compare with
1195 // sector = Int_t(track->GetAlpha() / AliTRDgeometry::GetAlpha());
1197 AliTRDseedV1 *tracklet = NULL;
1199 for(Int_t itl = 0; itl < AliTRDgeometry::kNlayer; itl++){
1200 if(!(tracklet = fkTrack->GetTracklet(itl)) || !tracklet->IsOK()) continue;
1201 sector = static_cast<Int_t>(tracklet->GetDetector()/AliTRDgeometry::kNdets);
1209 //________________________________________________________
1210 void AliTRDcheckDET::GetDistanceToTracklet(Double_t *dist, AliTRDseedV1 * const tracklet, AliTRDcluster * const c)
1212 Float_t x = c->GetX();
1213 dist[0] = c->GetY() - tracklet->GetYat(x);
1214 dist[1] = c->GetZ() - tracklet->GetZat(x);
1217 //________________________________________________________
1218 void AliTRDcheckDET::GetEtaPhiAt(const AliExternalTrackParam *track, Double_t x, Double_t &eta, Double_t &phi){
1220 // Get phi and eta at a given radial position
1222 AliExternalTrackParam workpar(*track);
1224 Double_t posLocal[3];
1225 Bool_t sucPos = workpar.GetXYZAt(x, fEvent->GetRunInfo()->GetMagneticField(), posLocal);
1226 Double_t sagPhi = sucPos ? TMath::ATan2(posLocal[1], posLocal[0]) : 0.;
1228 eta = workpar.Eta();
1232 //_______________________________________________________
1233 TH1* AliTRDcheckDET::MakePlotChi2() const
1235 // Plot chi2/track normalized to number of degree of freedom
1236 // (tracklets) and compare with the theoretical distribution.
1238 // Alex Bercuci <A.Bercuci@gsi.de>
1242 /* TH2S *h2 = (TH2S*)fContainer->At(kChi2);
1243 TF1 f("fChi2", "[0]*pow(x, [1]-1)*exp(-0.5*x)", 0., 50.);
1244 f.SetParLimits(1,1, 1e100);
1245 TLegend *leg = new TLegend(.7,.7,.95,.95);
1246 leg->SetBorderSize(1); leg->SetHeader("Tracklets per Track");
1248 Bool_t kFIRST = kTRUE;
1249 for(Int_t il=1; il<=h2->GetNbinsX(); il++){
1250 h1 = h2->ProjectionY(Form("pyChi2%d", il), il, il);
1251 if(h1->Integral()<50) continue;
1252 h1->Scale(1./h1->Integral());
1253 h1->SetMarkerStyle(7);h1->SetMarkerColor(il);
1254 h1->SetLineColor(il);h1->SetLineStyle(2);
1255 f.SetParameter(1, .5*il);f.SetLineColor(il);
1256 h1->Fit(&f, "QW+", kFIRST ? "pc": "pcsame");
1257 leg->AddEntry(h1, Form("%d", il), "l");
1259 h1->GetXaxis()->SetRangeUser(0., 25.);
1269 //________________________________________________________
1270 TH1* AliTRDcheckDET::MakePlotNTracklets(){
1272 // Make nice bar plot of the number of tracklets in each method
1274 TH2F *tmp = (TH2F *)fContainer->FindObject("hNtlsBAR");
1275 TH1D *hBAR = tmp->ProjectionY();
1276 TH1F *hSTA = (TH1F *)fContainer->FindObject("hNtlsSTA");
1277 TH1F *hCON = (TH1F *)fContainer->FindObject("hNtls");
1278 TLegend *leg = new TLegend(0.13, 0.75, 0.39, 0.89);
1279 leg->SetBorderSize(1);
1280 leg->SetFillColor(0);
1282 Float_t scale = hCON->Integral();
1283 if(scale) hCON->Scale(100./scale);
1284 hCON->SetFillColor(kRed);hCON->SetLineColor(kRed);
1285 hCON->SetBarWidth(0.2);
1286 hCON->SetBarOffset(0.6);
1287 hCON->SetStats(kFALSE);
1288 hCON->GetYaxis()->SetRangeUser(0.,40.);
1289 hCON->GetYaxis()->SetTitleOffset(1.2);
1290 hCON->Draw("bar1"); leg->AddEntry(hCON, "Total", "f");
1291 hCON->SetMaximum(55.);
1293 if(scale) hBAR->Scale(100./scale);
1294 hBAR->SetFillColor(kGreen);hBAR->SetLineColor(kGreen);
1295 hBAR->SetBarWidth(0.2);
1296 hBAR->SetBarOffset(0.2);
1298 hBAR->SetStats(kFALSE);
1299 hBAR->GetYaxis()->SetRangeUser(0.,40.);
1300 hBAR->GetYaxis()->SetTitleOffset(1.2);
1301 hBAR->Draw("bar1same"); leg->AddEntry(hBAR, "Barrel", "f");
1303 if(scale) hSTA->Scale(100./scale);
1304 hSTA->SetFillColor(kBlue);hSTA->SetLineColor(kBlue);
1305 hSTA->SetBarWidth(0.2);
1306 hSTA->SetBarOffset(0.4);
1308 hSTA->SetStats(kFALSE);
1309 hSTA->GetYaxis()->SetRangeUser(0.,40.);
1310 hSTA->GetYaxis()->SetTitleOffset(1.2);
1311 hSTA->Draw("bar1same"); leg->AddEntry(hSTA, "Stand Alone", "f");
1317 //________________________________________________________
1318 void AliTRDcheckDET::MakePlotnTrackletsVsP(){
1320 // Plot abundance of tracks with number of tracklets as function of momentum
1326 Color_t color[AliTRDgeometry::kNlayer] = {kBlue, kOrange, kBlack, kGreen, kCyan, kRed};
1327 TH1 *h(NULL); TGraphErrors *g[AliTRDgeometry::kNlayer];
1328 for(Int_t itl(0); itl<AliTRDgeometry::kNlayer; itl++){
1329 g[itl] = new TGraphErrors();
1330 g[itl]->SetLineColor(color[itl]);
1331 g[itl]->SetMarkerColor(color[itl]);
1332 g[itl]->SetMarkerStyle(20 + itl);
1335 TH2 *hBar = (TH2F *)fContainer->FindObject("hNtlsBAR");
1336 TAxis *ax(hBar->GetXaxis());
1337 Int_t np(ax->GetNbins());
1338 for(Int_t ipBin(1); ipBin<np; ipBin++){
1339 h = hBar->ProjectionY("npBin", ipBin, ipBin);
1340 if(!Int_t(h->Integral())) continue;
1341 h->Scale(100./h->Integral());
1342 Float_t p(ax->GetBinCenter(ipBin));
1343 Float_t dp(ax->GetBinWidth(ipBin));
1344 Int_t ip(g[0]->GetN());
1345 for(Int_t itl(AliTRDgeometry::kNlayer); itl--;){
1346 g[itl]->SetPoint(ip, p, h->GetBinContent(itl+1));
1347 g[itl]->SetPointError(ip, dp/3.46, h->GetBinError(itl+1));
1351 TLegend *leg = new TLegend(0.76, 0.6, 1., 0.9);
1352 leg->SetBorderSize(0);
1353 leg->SetHeader("Tracklet/Track");
1354 leg->SetFillStyle(0);
1355 h = hBar->ProjectionX("npxBin"); h->Reset();
1357 h->GetYaxis()->SetRangeUser(1., 99.);
1358 h->GetYaxis()->SetMoreLogLabels();
1359 h->GetYaxis()->CenterTitle();
1360 h->GetYaxis()->SetTitleOffset(1.2);
1361 h->SetYTitle("Prob. [%]");
1362 h->GetXaxis()->SetRangeUser(0.4, 12.);
1363 h->GetXaxis()->SetMoreLogLabels();
1364 h->GetXaxis()->CenterTitle();
1366 for(Int_t itl(AliTRDgeometry::kNlayer); itl--;){
1368 leg->AddEntry(g[itl], Form("n = %d", itl+1),"pl");
1372 gPad->SetLogx();gPad->SetLogy();
1375 //________________________________________________________
1376 Bool_t AliTRDcheckDET::MakePlotPulseHeight(){
1378 // Create Plot of the Pluse Height Spectrum
1380 TCanvas *output = gPad->GetCanvas();
1384 TObjArray *arr = (TObjArray*)fContainer->FindObject("<PH>");
1385 h = (TH1F*)arr->At(0);
1386 if((Int_t)h->GetEntries()){
1387 h->SetMarkerStyle(24);
1388 h->SetMarkerColor(kBlack);
1389 h->SetLineColor(kBlack);
1390 h->GetYaxis()->SetTitleOffset(1.5);
1392 // Trending for the pulse height: plateau value, slope and timebin of the maximum
1393 TLinearFitter fit(1,"pol1");
1395 for(Int_t itime = 10; itime <= 20; itime++){
1396 time = Double_t(itime);
1397 Double_t err(h->GetBinError(itime + 1));
1398 if(err>1.e-10) fit.AddPoint(&time, h->GetBinContent(itime + 1), err);
1401 Double_t plateau = fit.GetParameter(0) + 12 * fit.GetParameter(1);
1402 Double_t slope = fit.GetParameter(1);
1403 PutTrendValue("PHplateau", plateau);
1404 PutTrendValue("PHslope", slope);
1405 PutTrendValue("PHamplificationPeak", static_cast<Double_t>(h->GetMaximumBin()-1));
1406 AliDebug(1, Form("plateau %f, slope %f, MaxTime %f", plateau, slope, static_cast<Double_t>(h->GetMaximumBin()-1)));
1409 // copy the second histogram in a new one with the same x-dimension as the phs with respect to time
1410 h1 = (TH1F *)arr->At(1);
1411 h2 = new TH1F("hphs1","Average PH", 31, -0.5, 30.5);
1412 for(Int_t ibin = h1->GetXaxis()->GetFirst(); ibin < h1->GetNbinsX(); ibin++)
1413 h2->SetBinContent(ibin, h1->GetBinContent(ibin));
1414 h2->SetMarkerStyle(22);
1415 h2->SetMarkerColor(kBlue);
1416 h2->SetLineColor(kBlue);
1420 // create axis according to the histogram dimensions of the original second histogram
1421 TGaxis *axis = new TGaxis(gPad->GetUxmin(),
1425 -0.08, 4.88, 510,"-L");
1426 axis->SetLineColor(kBlue);
1427 axis->SetLabelColor(kBlue);
1428 axis->SetTextColor(kBlue);
1429 axis->SetTitle("x_{0}-x_{c} [cm]");
1433 TH2 *ph2d = (TH2F *)arr->At(2);
1434 ph2d->GetYaxis()->SetTitleOffset(1.8);
1435 ph2d->SetStats(kFALSE);
1440 //________________________________________________________
1441 void AliTRDcheckDET::MakePlotMeanClustersLayer(){
1443 // Create Summary plot for the mean number of clusters per layer
1445 TCanvas *output = gPad->GetCanvas();
1446 output->Divide(3,2);
1447 TObjArray *histos = (TObjArray *)fContainer->At(kNclustersLayer);
1449 AliWarning("Histos for each layer not found");
1452 TProfile2D *hlayer = NULL;
1453 for(Int_t ily = 0; ily < AliTRDgeometry::kNlayer; ily++){
1454 if(!(hlayer = dynamic_cast<TProfile2D *>(histos->At(ily)))) continue;
1455 output->cd(ily + 1);
1457 hlayer->Draw("colz");
1461 //________________________________________________________
1462 Bool_t AliTRDcheckDET::MakeBarPlot(TH1 *histo, Int_t color){
1464 // Draw nice bar plots
1466 if(!histo->GetEntries()) return kFALSE;
1467 histo->Scale(100./histo->Integral());
1468 histo->SetFillColor(color);
1469 histo->SetBarOffset(.2);
1470 histo->SetBarWidth(.6);
1471 histo->Draw("bar1");