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