]> git.uio.no Git - u/mrichter/AliRoot.git/blame - TRD/qaRec/AliTRDcheckDetector.cxx
Fixes to the Kalman filter (Markus)
[u/mrichter/AliRoot.git] / TRD / qaRec / AliTRDcheckDetector.cxx
CommitLineData
95cda7cf 1#include <TFile.h>
23abf4db 2#include <TH1F.h>
abc70aaf 3#include <TGraph.h>
23abf4db 4#include <TMath.h>
2a4a428a 5#include <TMap.h>
23abf4db 6#include <TObjArray.h>
2a4a428a 7#include <TObject.h>
8#include <TObjString.h>
23abf4db 9#include <TProfile.h>
9e548ae7 10#include <TProfile2D.h>
95cda7cf 11#include <TROOT.h>
23abf4db 12
107fde80 13#include "AliLog.h"
23abf4db 14#include "AliTRDcluster.h"
a24151d1 15#include "AliESDHeader.h"
16#include "AliESDRun.h"
3cfaffa4 17#include "AliESDtrack.h"
2374fb85 18#include "AliTRDgeometry.h"
3cfaffa4 19#include "AliTRDpadPlane.h"
23abf4db 20#include "AliTRDseedV1.h"
21#include "AliTRDtrackV1.h"
3cfaffa4 22#include "AliTRDtrackerV1.h"
23#include "AliTRDReconstructor.h"
4e488f26 24#include "AliTrackReference.h"
3cfaffa4 25#include "AliTrackPointArray.h"
26#include "AliTracker.h"
23abf4db 27#include "TTreeStream.h"
28
29#include "AliTRDtrackInfo/AliTRDtrackInfo.h"
a24151d1 30#include "AliTRDtrackInfo/AliTRDeventInfo.h"
23abf4db 31#include "AliTRDcheckDetector.h"
32
95cda7cf 33#include <cstdio>
56f8e2b1 34#include <iostream>
95cda7cf 35
23abf4db 36////////////////////////////////////////////////////////////////////////////
37// //
38// Reconstruction QA //
39// //
40// Task doing basic checks for tracking and detector performance //
41// //
42// Authors: //
43// Anton Andronic <A.Andronic@gsi.de> //
44// Markus Fasel <M.Fasel@gsi.de> //
45// //
46////////////////////////////////////////////////////////////////////////////
47
48//_______________________________________________________
49AliTRDcheckDetector::AliTRDcheckDetector():
93e41bce 50 AliTRDrecoTask("DetChecker", "Basic Detector Checker")
a24151d1 51 ,fEventInfo(0x0)
2a4a428a 52 ,fTriggerNames(0x0)
3cfaffa4 53 ,fReconstructor(0x0)
23abf4db 54{
2b468513 55 //
56 // Default constructor
57 //
a24151d1 58 DefineInput(1,AliTRDeventInfo::Class());
3cfaffa4 59 fReconstructor = new AliTRDReconstructor;
60 fReconstructor->SetRecoParam(AliTRDrecoParam::GetLowFluxParam());
107fde80 61 InitFunctorList();
23abf4db 62}
63
64//_______________________________________________________
65AliTRDcheckDetector::~AliTRDcheckDetector(){
2b468513 66 //
67 // Destructor
68 //
a24151d1 69 if(fEventInfo) delete fEventInfo;
2a4a428a 70 if(fTriggerNames) delete fTriggerNames;
3cfaffa4 71 delete fReconstructor;
23abf4db 72}
73
a24151d1 74//_______________________________________________________
75void AliTRDcheckDetector::ConnectInputData(Option_t *opt){
a391a274 76 //
77 // Connect the Input data with the task
78 //
79 AliTRDrecoTask::ConnectInputData(opt);
80 fEventInfo = dynamic_cast<AliTRDeventInfo *>(GetInputData(1));
a24151d1 81}
82
23abf4db 83//_______________________________________________________
84void AliTRDcheckDetector::CreateOutputObjects(){
2b468513 85 //
86 // Create Output Objects
87 //
2374fb85 88 OpenFile(0,"RECREATE");
107fde80 89 fContainer = Histos();
2374fb85 90 if(!fTriggerNames) fTriggerNames = new TMap();
23abf4db 91}
92
93//_______________________________________________________
107fde80 94void AliTRDcheckDetector::Exec(Option_t *opt){
2b468513 95 //
96 // Execution function
97 // Filling TRD quality histos
98 //
9e548ae7 99 if(!HasMCdata() && fEventInfo->GetEventHeader()->GetEventType() != 7) return; // For real data we select only physical events
107fde80 100 AliTRDrecoTask::Exec(opt);
2b468513 101 Int_t nTracks = 0; // Count the number of tracks per event
a24151d1 102 Int_t triggermask = fEventInfo->GetEventHeader()->GetTriggerMask();
afefec95 103 TString triggername = fEventInfo->GetRunInfo()->GetFiredTriggerClasses(triggermask);
9e548ae7 104 if(fDebugLevel > 6)printf("Trigger cluster: %d, Trigger class: %s\n", triggermask, triggername.Data());
a24151d1 105 dynamic_cast<TH1F *>(fContainer->UncheckedAt(kNEventsTrigger))->Fill(triggermask);
2b468513 106 for(Int_t iti = 0; iti < fTracks->GetEntriesFast(); iti++){
107fde80 107 if(!fTracks->UncheckedAt(iti)) continue;
108 AliTRDtrackInfo *fTrackInfo = dynamic_cast<AliTRDtrackInfo *>(fTracks->UncheckedAt(iti));
109 if(!fTrackInfo->GetTrack()) continue;
2b468513 110 nTracks++;
111 }
959f550d 112 if(nTracks){
a391a274 113 dynamic_cast<TH1F *>(fContainer->UncheckedAt(kNEventsTriggerTracks))->Fill(triggermask);
114 dynamic_cast<TH1F *>(fContainer->UncheckedAt(kNTracksEventHist))->Fill(nTracks);
959f550d 115 }
9e548ae7 116 if(triggermask <= 20 && !fTriggerNames->FindObject(Form("%d", triggermask))){
2a4a428a 117 fTriggerNames->Add(new TObjString(Form("%d", triggermask)), new TObjString(triggername));
9e548ae7 118 // also set the label for both histograms
119 TH1 *histo = dynamic_cast<TH1F *>(fContainer->UncheckedAt(kNEventsTriggerTracks));
120 histo->GetXaxis()->SetBinLabel(histo->FindBin(triggermask), triggername);
121 histo = dynamic_cast<TH1F *>(fContainer->UncheckedAt(kNEventsTrigger));
122 histo->GetXaxis()->SetBinLabel(histo->FindBin(triggermask), triggername);
123 }
2b468513 124 PostData(0, fContainer);
23abf4db 125}
126
127//_______________________________________________________
128void AliTRDcheckDetector::Terminate(Option_t *){
2b468513 129 //
130 // Terminate function
131 //
23abf4db 132}
133
95cda7cf 134//_______________________________________________________
135Bool_t AliTRDcheckDetector::PostProcess(){
a391a274 136 //
137 // Do Postprocessing (for the moment set the number of Reference histograms)
138 //
139
140 TH1 * histo = 0x0;
141 histo = dynamic_cast<TH1F *>(fContainer->UncheckedAt(kNTracksEventHist));
142 histo->GetXaxis()->SetTitle("Number of Tracks");
143 histo->GetYaxis()->SetTitle("Events");
144 histo = dynamic_cast<TH1F *>(fContainer->UncheckedAt(kNclustersHist));
145 histo->GetXaxis()->SetTitle("Number of Clusters");
146 histo->GetYaxis()->SetTitle("Entries");
147 histo = dynamic_cast<TH1F *>(fContainer->UncheckedAt(kNtrackletsHist));
148 histo->GetXaxis()->SetTitle("Number of Tracklets");
149 histo->GetYaxis()->SetTitle("Entries");
150 histo = dynamic_cast<TH1F *>(fContainer->UncheckedAt(kNclusterTrackletHist));
151 histo->GetXaxis()->SetTitle("Number of Clusters");
152 histo->GetYaxis()->SetTitle("Entries");
153 histo = dynamic_cast<TH1F *>(fContainer->UncheckedAt(kChi2));
154 histo->GetXaxis()->SetTitle("#chi^2");
155 histo->GetYaxis()->SetTitle("Entries");
156 histo = dynamic_cast<TH1F *>(fContainer->UncheckedAt(kNTracksSectorHist));
157 histo->GetXaxis()->SetTitle("Sector");
158 histo->GetYaxis()->SetTitle("Number of Tracks");
159 histo = dynamic_cast<TProfile *>(fContainer->UncheckedAt(kPulseHeight));
160 histo->GetXaxis()->SetTitle("Time / 100ns");
161 histo->GetYaxis()->SetTitle("Average Pulse Height (a. u.)");
162 histo = dynamic_cast<TH1F *>(fContainer->UncheckedAt(kClusterCharge));
163 histo->GetXaxis()->SetTitle("Cluster Charge (a.u.)");
164 histo->GetYaxis()->SetTitle("Entries");
165 histo = dynamic_cast<TH1F *>(fContainer->UncheckedAt(kChargeDeposit));
166 histo->GetXaxis()->SetTitle("Charge Deposit (a.u.)");
167 histo->GetYaxis()->SetTitle("Entries");
168
169 // Calculate the purity of the trigger clusters
170 histo = dynamic_cast<TH1F *>(fContainer->UncheckedAt(kNEventsTrigger));
171 TH1F *histoTracks = dynamic_cast<TH1F *>(fContainer->UncheckedAt(kNEventsTriggerTracks));
172 histoTracks->Divide(histo);
173 Float_t purities[20], val = 0;
174 TString triggernames[20];
175 Int_t nTriggerClasses = 0;
176 for(Int_t ibin = 1; ibin <= histo->GetNbinsX(); ibin++){
177 if((val = histoTracks->GetBinContent(ibin))){
178 purities[nTriggerClasses] = val;
179 triggernames[nTriggerClasses] = histoTracks->GetXaxis()->GetBinLabel(ibin);
180 nTriggerClasses++;
181 }
182 }
183 TH1F *hTriggerInf = new TH1F("fTriggerInf", "Trigger Information", TMath::Max(nTriggerClasses, 1), 0, TMath::Max(nTriggerClasses, 1));
184 for(Int_t ibin = 1; ibin <= nTriggerClasses; ibin++){
185 hTriggerInf->SetBinContent(ibin, purities[ibin-1]);
186 hTriggerInf->GetXaxis()->SetBinLabel(ibin, triggernames[ibin-1].Data());
187 }
188 hTriggerInf->GetXaxis()->SetTitle("Trigger Cluster");
189 hTriggerInf->GetYaxis()->SetTitle("Ratio");
190 hTriggerInf->GetYaxis()->SetRangeUser(0,1);
9e548ae7 191// hTriggerInf->SetMarkerColor(kBlue);
192// hTriggerInf->SetMarkerStyle(22);
a391a274 193 fContainer->Add(hTriggerInf);
194 fNRefFigures = 10;
195 return kTRUE;
95cda7cf 196}
197
198//_______________________________________________________
a391a274 199void AliTRDcheckDetector::GetRefFigure(Int_t ifig){
200 //
201 // Setting Reference Figures
202 //
203 TH1 *h = 0x0;
204 switch(ifig){
205 case 0:
206 ((TH1F*)fContainer->At(kNTracksEventHist))->Draw("pl");
207 break;
208 case 1:
209 ((TH1F*)fContainer->At(kNclustersHist))->Draw("pl");
210 break;
211 case 2:
212 h = (TH1F*)fContainer->At(kNtrackletsHist);
213 if(!h->GetEntries()) break;
214 h->Scale(100./h->Integral());
215 h->GetXaxis()->SetRangeUser(.5, 6.5);
216 h->SetFillColor(kGreen);
217 h->SetBarOffset(.2);
218 h->SetBarWidth(.6);
219 h->Draw("bar1");
220 break;
221 case 3:
222 ((TH1F*)fContainer->At(kNclusterTrackletHist))->Draw("pc");
223 break;
224 case 4:
225 ((TH1F*)fContainer->At(kChi2))->Draw("");
226 break;
227 case 5:
228 h = (TH1F*)fContainer->At(kNTracksSectorHist);
229 if(!h->GetEntries()) break;
230 h->Scale(100./h->Integral());
231 h->SetFillColor(kGreen);
232 h->SetBarOffset(.2);
233 h->SetBarWidth(.6);
234 h->Draw("bar1");
235 break;
236 case 6:
237 h = (TH1F*)fContainer->At(kPulseHeight);
238 h->SetMarkerStyle(24);
239 h->Draw("e1");
240 break;
241 case 7:
242 ((TH1F*)fContainer->At(kClusterCharge))->Draw("c");
243 break;
244 case 8:
245 ((TH1F*)fContainer->At(kChargeDeposit))->Draw("c");
246 break;
247 case 9:
248 h=(TH1F*)fContainer->At(kPurity);
249 h->SetBarOffset(.2);
250 h->SetBarWidth(.6);
251 h->SetFillColor(kGreen);
252 h->Draw("bar1");
253 break;
254 default:
255 ((TH1F*)fContainer->At(kNTracksEventHist))->Draw("pl");
256 break;
257 }
95cda7cf 258}
259
107fde80 260//_______________________________________________________
261TObjArray *AliTRDcheckDetector::Histos(){
262 //
263 // Create QA histograms
264 //
265 if(fContainer) return fContainer;
266
267 fContainer = new TObjArray(25);
268 // Register Histograms
269 fContainer->AddAt(new TH1F("hNtrks", "Number of Tracks per event", 100, 0, 100), kNTracksEventHist);
270 fContainer->AddAt(new TH1F("hEventsTriggerTracks", "Trigger Class (Tracks)", 100, 0, 100), kNEventsTriggerTracks);
271 fContainer->AddAt(new TH1F("hNcls", "Nr. of clusters per track", 181, -0.5, 180.5), kNclustersHist);
272 fContainer->AddAt(new TH1F("hNtls", "Nr. tracklets per track", 7, -0.5, 6.5), kNtrackletsHist);
3cfaffa4 273 fContainer->AddAt(new TH1F("hNtlsFindable", "Ratio of found/findable Tracklets" , 11, -0.05, 1.05), kNTrackletsVsFindable);
107fde80 274 fContainer->AddAt(new TH1F("hNclTls","Mean Number of clusters per tracklet", 31, -0.5, 30.5), kNclusterTrackletHist);
275 fContainer->AddAt(new TH1F("hChi2", "Chi2", 200, 0, 20), kChi2);
276 fContainer->AddAt(new TH1F("hChi2n", "Norm. Chi2 (tracklets)", 50, 0, 5), kChi2Normalized);
277 fContainer->AddAt(new TH1F("hSM", "Track Counts in Supermodule", 18, -0.5, 17.5), kNTracksSectorHist);
278 // Detector signal on Detector-by-Detector basis
279 fContainer->AddAt(new TProfile("hPHdetector", "Average PH", 31, -0.5, 30.5), kPulseHeight);
280 fContainer->AddAt(new TH1F("hQclDetector", "Cluster charge", 200, 0, 1200), kClusterCharge);
281 fContainer->AddAt(new TH1F("hQTdetector", "Total Charge Deposit", 6000, 0, 6000), kChargeDeposit);
282 fContainer->AddAt(new TH1F("hEventsTrigger", "Trigger Class", 100, 0, 100), kNEventsTrigger);
283
284 return fContainer;
285}
286
287/*
288* Plotting Functions
289*/
290
291//_______________________________________________________
2374fb85 292TH1 *AliTRDcheckDetector::PlotMeanNClusters(const AliTRDtrackV1 *track){
107fde80 293 //
294 // Plot the mean number of clusters per tracklet
295 //
74b2e03d 296 if(track) fTrack = track;
107fde80 297 if(!fTrack){
74b2e03d 298 AliWarning("No Track defined.");
299 return 0x0;
107fde80 300 }
301 TH1 *h = 0x0;
302 if(!(h = dynamic_cast<TH1F *>(fContainer->At(kNclusterTrackletHist)))){
303 AliWarning("No Histogram defined.");
304 return 0x0;
305 }
306 AliTRDseedV1 *tracklet = 0x0;
2374fb85 307 for(Int_t itl = 0; itl < AliTRDgeometry::kNlayer; itl++){
107fde80 308 if(!(tracklet = fTrack->GetTracklet(itl)) || !tracklet->IsOK()) continue;
309 h->Fill(tracklet->GetN());
310 }
311 return h;
312}
313
314//_______________________________________________________
2374fb85 315TH1 *AliTRDcheckDetector::PlotNClusters(const AliTRDtrackV1 *track){
107fde80 316 //
317 // Plot the number of clusters in one track
318 //
74b2e03d 319 if(track) fTrack = track;
107fde80 320 if(!fTrack){
74b2e03d 321 AliWarning("No Track defined.");
322 return 0x0;
107fde80 323 }
324 TH1 *h = 0x0;
325 if(!(h = dynamic_cast<TH1F *>(fContainer->At(kNclustersHist)))){
326 AliWarning("No Histogram defined.");
327 return 0x0;
328 }
2374fb85 329
107fde80 330 Int_t nclusters = 0;
331 AliTRDseedV1 *tracklet = 0x0;
2374fb85 332 for(Int_t itl = 0; itl < AliTRDgeometry::kNlayer; itl++){
107fde80 333 if(!(tracklet = fTrack->GetTracklet(itl)) || !tracklet->IsOK()) continue;
334 nclusters += tracklet->GetN();
335 if(fDebugLevel > 2){
336 Int_t crossing = tracklet->GetNChange();
337 AliTRDcluster *c = 0x0;
3cfaffa4 338 for(Int_t itime = 0; itime < AliTRDtrackerV1::GetNTimeBins(); itime++){
107fde80 339 if(!(c = tracklet->GetClusters(itime))) continue;
340 break;
341 }
342 Int_t detector = c->GetDetector();
2374fb85 343 Float_t sector = static_cast<Int_t>(detector/AliTRDgeometry::kNdets);
107fde80 344 Float_t theta = TMath::ATan(tracklet->GetZfit(1));
345 Float_t phi = TMath::ATan(tracklet->GetYfit(1));
346 Float_t momentum = 0.;
347 Int_t pdg = 0;
52a79f8d 348 Int_t kinkIndex = fESD ? fESD->GetKinkIndex() : 0;
349 UShort_t TPCncls = fESD ? fESD->GetTPCncls() : 0;
107fde80 350 if(fMC){
52a79f8d 351 if(fMC->GetTrackRef()) momentum = fMC->GetTrackRef()->P();
75d00a6d 352 pdg = fMC->GetPDG();
107fde80 353 }
354 (*fDebugStream) << "NClusters"
355 << "Detector=" << detector
356 << "Sector=" << sector
357 << "crossing=" << crossing
358 << "momentum=" << momentum
359 << "pdg=" << pdg
360 << "theta=" << theta
361 << "phi=" << phi
52a79f8d 362 << "kinkIndex=" << kinkIndex
363 << "TPCncls=" << TPCncls
107fde80 364 << "nclusters=" << nclusters
365 << "\n";
366 }
367 }
368 h->Fill(nclusters);
369 return h;
370}
371
372//_______________________________________________________
2374fb85 373TH1 *AliTRDcheckDetector::PlotChi2(const AliTRDtrackV1 *track){
107fde80 374 //
375 // Plot the chi2 of the track
376 //
74b2e03d 377 if(track) fTrack = track;
378 if(!fTrack){
379 AliWarning("No Track defined.");
380 return 0x0;
107fde80 381 }
382 TH1 *h = 0x0;
383 if(!(h = dynamic_cast<TH1F *>(fContainer->At(kChi2)))){
384 AliWarning("No Histogram defined.");
385 return 0x0;
386 }
387 h->Fill(fTrack->GetChi2());
388 return h;
389}
390
391//_______________________________________________________
2374fb85 392TH1 *AliTRDcheckDetector::PlotNormalizedChi2(const AliTRDtrackV1 *track){
107fde80 393 //
394 // Plot the chi2 of the track
395 //
74b2e03d 396 if(track) fTrack = track;
397 if(!fTrack){
398 AliWarning("No Track defined.");
399 return 0x0;
107fde80 400 }
401 TH1 *h = 0x0;
402 if(!(h = dynamic_cast<TH1F *>(fContainer->At(kChi2Normalized)))){
403 AliWarning("No Histogram defined.");
404 return 0x0;
405 }
107fde80 406 Int_t nTracklets = 0;
2374fb85 407 AliTRDseedV1 *tracklet = 0x0;
408 for(Int_t itl = 0; itl < AliTRDgeometry::kNlayer; itl++){
409 if(!(tracklet = fTrack->GetTracklet(itl)) || !tracklet->IsOK()) continue;
410 nTracklets++;
411 }
107fde80 412 h->Fill(fTrack->GetChi2()/nTracklets);
413 return h;
414}
415
416
417//_______________________________________________________
2374fb85 418TH1 *AliTRDcheckDetector::PlotNTracklets(const AliTRDtrackV1 *track){
107fde80 419 //
420 // Plot the number of tracklets
421 //
74b2e03d 422 if(track) fTrack = track;
107fde80 423 if(!fTrack){
74b2e03d 424 AliWarning("No Track defined.");
425 return 0x0;
107fde80 426 }
427 TH1 *h = 0x0;
428 if(!(h = dynamic_cast<TH1F *>(fContainer->At(kNtrackletsHist)))){
429 AliWarning("No Histogram defined.");
430 return 0x0;
431 }
3cfaffa4 432 Int_t nTracklets = GetNTracklets(track);
433 h->Fill(nTracklets);
434 if(fDebugLevel > 3){
435 if(nTracklets == 1){
436 // If we have one Tracklet, check in which layer this happens
437 Int_t layer = -1;
438 AliTRDseedV1 *tracklet = 0x0;
439 for(Int_t il = 0; il < AliTRDgeometry::kNlayer; il++){
440 if((tracklet = track->GetTracklet(il)) && tracklet->IsOK()){layer = il; break;}
441 }
442 (*fDebugStream) << "PlotNTracklets"
443 << "Layer=" << layer
444 << "\n";
445 }
446 }
447 return h;
448}
449
450//_______________________________________________________
451TH1 *AliTRDcheckDetector::PlotTrackletsVsFindable(const AliTRDtrackV1 *track){
452 //
453 // Plots the ratio of number of tracklets vs.
454 // number of findable tracklets
455 //
456 // Findable tracklets are defined as track prolongation
457 // to layer i does not hit the dead area +- epsilon
458 //
459 const Float_t epsilon = 0.01; // dead area tolerance
460 const Float_t epsilon_R = 1;
461 Double_t x_anode[AliTRDgeometry::kNlayer] = {300.2, 312.8, 325.4, 338.0, 350.6, 363.2}; // Take the default X0
462
463 if(track) fTrack = track;
464 if(!fTrack){
465 AliWarning("No Track defined.");
466 return 0x0;
467 }
468 TH1 *h = 0x0;
469 if(!(h = dynamic_cast<TH1F *>(fContainer->At(kNTrackletsVsFindable)))){
470 AliWarning("No Histogram defined.");
471 return 0x0;
472 }
473 Int_t nFound = 0, nFindable = 0;
474 Int_t stack = -1;
475 Double_t ymin = 0., ymax = 0., zmin = 0., zmax = 0.;
476 Double_t y = 0., z = 0.;
107fde80 477 AliTRDseedV1 *tracklet = 0x0;
3cfaffa4 478 AliTRDpadPlane *pp;
479 AliTRDgeometry *fGeo = new AliTRDgeometry;
480 for(Int_t il = 0; il < AliTRDgeometry::kNlayer; il++){
481 if((tracklet = fTrack->GetTracklet(il)) && tracklet->IsOK()){
482 tracklet->SetReconstructor(fReconstructor);
483 nFound++;
484 // printout tracklet
485 }
107fde80 486 }
3cfaffa4 487 // 2 Different cases:
488 // 1st stand alone: here we cannot propagate, but be can do a Tilted Rieman Fit
489 // 2nd barrel track: here we propagate the track to the layers
490 AliTrackPoint points[6];
491 Float_t xyz[3];
492 memset(xyz, 0, sizeof(Float_t) * 3);
493 if(((fESD->GetStatus() & AliESDtrack::kTRDout) > 0) && !((fESD->GetStatus() & AliESDtrack::kTRDin) > 0)){
494 // stand alone track
495 for(Int_t il = 0; il < AliTRDgeometry::kNlayer; il++){
496 xyz[0] = x_anode[il];
497 points[il].SetXYZ(xyz);
498 }
499 AliTRDtrackerV1::FitRiemanTilt(const_cast<AliTRDtrackV1 *>(fTrack), 0x0, kTRUE, 6, points);
500 } else {
501 // barrel track
502 for(Int_t il = AliTRDgeometry::kNlayer; il--;){
503 // The track points have to be in reverse order for the Kalman Filter (since we move back)
504 xyz[0] = x_anode[AliTRDgeometry::kNlayer - il - 1] - epsilon_R;
505 points[il].SetXYZ(xyz);
506 }
507 for(Int_t ipt = 0; ipt < AliTRDgeometry::kNlayer; ipt++)
508 printf("%d. X = %f\n", ipt, points[ipt].GetX());
509 AliTRDtrackerV1::FitKalman(const_cast<AliTRDtrackV1 *>(fTrack), 0x0, kFALSE, 6, points);
510 AliTrackPoint tempPointCont[AliTRDgeometry::kNlayer];
511 memcpy(tempPointCont, points, sizeof(AliTrackPoint) * AliTRDgeometry::kNlayer);
512 for(Int_t il = AliTRDgeometry::kNlayer; il--; ){
513 tempPointCont[il].GetXYZ(xyz);
514 points[AliTRDgeometry::kNlayer - il -1].SetXYZ(xyz);
515 }
516 }
517 for(Int_t il = 0; il < AliTRDgeometry::kNlayer; il++){
518 y = points[il].GetY();
519 z = points[il].GetZ();
520 if((stack = fGeo->GetStack(z, il)) < 0) continue; // Not findable
521 pp = fGeo->GetPadPlane(il, stack);
522 ymin = pp->GetCol0() + epsilon;
523 ymax = pp->GetColEnd() - epsilon;
524 zmin = pp->GetRowEnd() + epsilon;
525 zmax = pp->GetRow0() - epsilon;
526 // ignore y-crossing (material)
527 if((z > zmin && z < zmax) && (y > ymin && y < ymax)) nFindable++;
528/*
529 if(fDebugLevel > 3){
530 Double_t pos_tracklet[2] = {tracklet ? tracklet->GetYfit(0) : 0, tracklet ? tracklet->GetMeanz() : 0};
531 Int_t hasTracklet = tracklet ? 1 : 0;
532 (*fDebugStream) << "GetFindableTracklets"
533 << "layer=" << il
534 << "ytracklet=" << pos_tracklet[0]
535 << "ytrack=" << y
536 << "ztracklet=" << pos_tracklet[1]
537 << "ztrack=" << z
538 << "tracklet=" << hasTracklet
539 << "\n";
540 }
541*/
542 }
543
544 delete fGeo;
545 h->Fill(nFindable > 0 ? TMath::Min(nFound/static_cast<Double_t>(nFindable), 1.) : 1);
546 if(fDebugLevel > 2) AliInfo(Form("Findable[Found]: %d[%d|%f]", nFindable, nFound, nFound/static_cast<Float_t>(nFindable > 0 ? nFindable : 1)));
107fde80 547 return h;
548}
549
550//_______________________________________________________
2374fb85 551TH1 *AliTRDcheckDetector::PlotPulseHeight(const AliTRDtrackV1 *track){
107fde80 552 //
553 // Plot the average pulse height
554 //
74b2e03d 555 if(track) fTrack = track;
107fde80 556 if(!fTrack){
74b2e03d 557 AliWarning("No Track defined.");
558 return 0x0;
107fde80 559 }
560 TProfile *h = 0x0;
561 if(!(h = dynamic_cast<TProfile *>(fContainer->At(kPulseHeight)))){
562 AliWarning("No Histogram defined.");
563 return 0x0;
564 }
107fde80 565 AliTRDseedV1 *tracklet = 0x0;
566 AliTRDcluster *c = 0x0;
2374fb85 567 for(Int_t itl = 0; itl < AliTRDgeometry::kNlayer; itl++){
107fde80 568 if(!(tracklet = fTrack->GetTracklet(itl)) || !tracklet->IsOK())continue;
3cfaffa4 569 for(Int_t itime = 0; itime < AliTRDtrackerV1::GetNTimeBins(); itime++){
107fde80 570 if(!(c = tracklet->GetClusters(itime))) continue;
571 Int_t localtime = c->GetLocalTimeBin();
572 Double_t absolute_charge = TMath::Abs(c->GetQ());
573 h->Fill(localtime, absolute_charge);
574 if(fDebugLevel > 3){
575 Int_t crossing = tracklet->GetNChange();
576 Int_t detector = c->GetDetector();
2374fb85 577 Float_t sector = static_cast<Int_t>(detector/AliTRDgeometry::kNdets);
107fde80 578 Float_t theta = TMath::ATan(tracklet->GetZfit(1));
579 Float_t phi = TMath::ATan(tracklet->GetYfit(1));
580 Float_t momentum = 0.;
581 Int_t pdg = 0;
52a79f8d 582 Int_t kinkIndex = fESD ? fESD->GetKinkIndex() : 0;
583 UShort_t TPCncls = fESD ? fESD->GetTPCncls() : 0;
107fde80 584 if(fMC){
75d00a6d 585 if(fMC->GetTrackRef()) momentum = fMC->GetTrackRef()->P();
586 pdg = fMC->GetPDG();
107fde80 587 }
588 (*fDebugStream) << "PulseHeight"
589 << "Detector=" << detector
590 << "Sector=" << sector
2374fb85 591 << "crossing=" << crossing
107fde80 592 << "Timebin=" << localtime
593 << "Charge=" << absolute_charge
594 << "momentum=" << momentum
595 << "pdg=" << pdg
596 << "theta=" << theta
597 << "phi=" << phi
52a79f8d 598 << "kinkIndex=" << kinkIndex
599 << "TPCncls=" << TPCncls
107fde80 600 << "\n";
601 }
602 }
603 }
604 return h;
605}
606
607//_______________________________________________________
2374fb85 608TH1 *AliTRDcheckDetector::PlotClusterCharge(const AliTRDtrackV1 *track){
107fde80 609 //
610 // Plot the cluster charge
611 //
74b2e03d 612 if(track) fTrack = track;
107fde80 613 if(!fTrack){
74b2e03d 614 AliWarning("No Track defined.");
615 return 0x0;
107fde80 616 }
617 TH1 *h = 0x0;
618 if(!(h = dynamic_cast<TH1F *>(fContainer->At(kClusterCharge)))){
619 AliWarning("No Histogram defined.");
620 return 0x0;
621 }
107fde80 622 AliTRDseedV1 *tracklet = 0x0;
623 AliTRDcluster *c = 0x0;
2374fb85 624 for(Int_t itl = 0; itl < AliTRDgeometry::kNlayer; itl++){
107fde80 625 if(!(tracklet = fTrack->GetTracklet(itl)) || !tracklet->IsOK())continue;
3cfaffa4 626 for(Int_t itime = 0; itime < AliTRDtrackerV1::GetNTimeBins(); itime++){
107fde80 627 if(!(c = tracklet->GetClusters(itime))) continue;
628 h->Fill(c->GetQ());
629 }
630 }
631 return h;
632}
633
634//_______________________________________________________
2374fb85 635TH1 *AliTRDcheckDetector::PlotChargeDeposit(const AliTRDtrackV1 *track){
107fde80 636 //
637 // Plot the charge deposit per chamber
638 //
74b2e03d 639 if(track) fTrack = track;
107fde80 640 if(!fTrack){
74b2e03d 641 AliWarning("No Track defined.");
642 return 0x0;
107fde80 643 }
644 TH1 *h = 0x0;
645 if(!(h = dynamic_cast<TH1F *>(fContainer->At(kChargeDeposit)))){
646 AliWarning("No Histogram defined.");
647 return 0x0;
648 }
107fde80 649 AliTRDseedV1 *tracklet = 0x0;
2374fb85 650 AliTRDcluster *c = 0x0, *c1 = 0x0; // c1 for the Debug Stream
107fde80 651 Double_t Qtot = 0;
3cfaffa4 652 Int_t nTracklets = 0;
653 if(fDebugLevel > 3)
654 nTracklets = GetNTracklets(fTrack); // fill NTracklet to the Debug Stream
2374fb85 655 for(Int_t itl = 0x0; itl < AliTRDgeometry::kNlayer; itl++){
107fde80 656 if(!(tracklet = fTrack->GetTracklet(itl)) || !tracklet->IsOK()) continue;
657 Qtot = 0;
2374fb85 658 c1 = 0x0;
3cfaffa4 659 for(Int_t itime = 0; itime < AliTRDtrackerV1::GetNTimeBins(); itime++){
107fde80 660 if(!(c = tracklet->GetClusters(itime))) continue;
2374fb85 661 if(!c1) c1 = c;
107fde80 662 Qtot += TMath::Abs(c->GetQ());
663 }
664 h->Fill(Qtot);
665 if(fDebugLevel > 3){
666 Int_t crossing = tracklet->GetNChange();
2374fb85 667 Int_t detector = c1->GetDetector();
668 Float_t sector = static_cast<Int_t>(detector/AliTRDgeometry::kNdets);
107fde80 669 Float_t theta = TMath::ATan(tracklet->GetZfit(1));
670 Float_t phi = TMath::ATan(tracklet->GetYfit(1));
671 Float_t momentum = 0.;
672 Int_t pdg = 0;
52a79f8d 673 Int_t kinkIndex = fESD ? fESD->GetKinkIndex() : 0;
674 UShort_t TPCncls = fESD ? fESD->GetTPCncls() : 0;
107fde80 675 if(fMC){
75d00a6d 676 if(fMC->GetTrackRef()) momentum = fMC->GetTrackRef()->P();
677 pdg = fMC->GetPDG();
107fde80 678 }
679 (*fDebugStream) << "ChargeDeposit"
680 << "Detector=" << detector
681 << "Sector=" << sector
682 << "crossing=" << crossing
683 << "momentum=" << momentum
3cfaffa4 684 << "nTracklets="<< nTracklets
107fde80 685 << "pdg=" << pdg
686 << "theta=" << theta
687 << "phi=" << phi
52a79f8d 688 << "kinkIndex=" << kinkIndex
689 << "TPCncls=" << TPCncls
107fde80 690 << "QT=" << Qtot
691 << "\n";
692 }
693 }
694 return h;
695}
696
697//_______________________________________________________
2374fb85 698TH1 *AliTRDcheckDetector::PlotTracksSector(const AliTRDtrackV1 *track){
107fde80 699 //
700 // Plot the number of tracks per Sector
701 //
74b2e03d 702 if(track) fTrack = track;
107fde80 703 if(!fTrack){
74b2e03d 704 AliWarning("No Track defined.");
705 return 0x0;
107fde80 706 }
707 TH1 *h = 0x0;
708 if(!(h = dynamic_cast<TH1F *>(fContainer->At(kNTracksSectorHist)))){
709 AliWarning("No Histogram defined.");
710 return 0x0;
711 }
107fde80 712 AliTRDseedV1 *tracklet = 0x0;
713 AliTRDcluster *c = 0x0;
714 Int_t sector = -1;
2374fb85 715 for(Int_t itl = 0; itl < AliTRDgeometry::kNlayer; itl++){
716 if(!(tracklet = fTrack->GetTracklet(itl)) || !tracklet->IsOK()) continue;
3cfaffa4 717 for(Int_t itime = 0; itime < AliTRDtrackerV1::GetNTimeBins(); itime++){
107fde80 718 if(!(c = tracklet->GetClusters(itime))) continue;
2374fb85 719 sector = static_cast<Int_t>(c->GetDetector()/AliTRDgeometry::kNdets);
107fde80 720 }
721 break;
722 }
723 h->Fill(sector);
724 return h;
725}
726
3cfaffa4 727//_______________________________________________________
728Int_t AliTRDcheckDetector::GetNTracklets(const AliTRDtrackV1 *track){
729 //
730 // Count the number of tracklets per track
731 //
732 if(!track) return 0;
733 Int_t nTracklets = 0;
734 AliTRDseedV1 *tracklet = 0x0;
735 for(Int_t il = 0; il < AliTRDgeometry::kNlayer; il++){
736 if((tracklet = track->GetTracklet(il)) && tracklet->IsOK()) nTracklets++;
737 }
738 return nTracklets;
739}
740
741//________________________________________________________
742void AliTRDcheckDetector::SetRecoParam(AliTRDrecoParam *r)
743{
744
745 fReconstructor->SetRecoParam(r);
746}