Coverity fix
[u/mrichter/AliRoot.git] / TRD / AliTRDQADataMaker.cxx
1 /**************************************************************************
2  * Copyright(c) 1998-1999, ALICE Experiment at CERN, All rights reserved. *
3  *                                                                        *
4  * Author: The ALICE Off-line Project.                                    *
5  * Contributors are mentioned in the code where appropriate.              *
6  *                                                                        *
7  * Permission to use, copy, modify and distribute this software and its   *
8  * documentation strictly for non-commercial purposes is hereby granted   *
9  * without fee, provided that the above copyright notice appears in all   *
10  * copies and that both the copyright notice and this permission notice   *
11  * appear in the supporting documentation. The authors make no claims     *
12  * about the suitability of this software for any purpose. It is          *
13  * provided "as is" without express or implied warranty.                  *
14  **************************************************************************/
15
16 /* $Id$ */
17
18 ////////////////////////////////////////////////////////////////////////////
19 //                                                                        //
20 //  Produces the data needed to calculate the quality assurance.          //
21 //  All data must be mergeable objects.                                   //
22 //                                                                        //
23 //  Author:                                                               //
24 //    Sylwester Radomski (radomski@physi.uni-heidelberg.de)               //
25 //                                                                        //
26 ////////////////////////////////////////////////////////////////////////////
27
28 // --- ROOT system ---
29 #include <TClonesArray.h>
30 #include <TH1F.h> 
31 #include <TH2F.h>
32 #include <TH3F.h>
33 #include <TProfile.h>
34 #include <TF1.h>
35
36 // --- AliRoot header files ---
37 #include "AliESDEvent.h"
38 #include "AliQAChecker.h"
39
40 #include "AliTRDdigit.h"
41 #include "AliTRDhit.h"
42 #include "AliTRDcluster.h"
43 #include "AliTRDQADataMaker.h"
44 #include "AliTRDdigitsManager.h"
45 #include "AliTRDgeometry.h"
46 #include "AliTRDarrayADC.h"
47
48 ClassImp(AliTRDQADataMaker)
49
50 //____________________________________________________________________________ 
51   AliTRDQADataMaker::AliTRDQADataMaker() : 
52   AliQADataMaker(AliQAv1::GetDetName(AliQAv1::kTRD), "TRD Quality Assurance Data Maker")
53 {
54   //
55   // Default constructor
56 }
57
58 //____________________________________________________________________________ 
59 AliTRDQADataMaker::AliTRDQADataMaker(const AliTRDQADataMaker& qadm) :
60   AliQADataMaker()
61 {
62   //
63   // Copy constructor 
64   //
65
66   SetName((const char*)qadm.GetName()) ; 
67   SetTitle((const char*)qadm.GetTitle()); 
68
69 }
70
71 //__________________________________________________________________
72 AliTRDQADataMaker& AliTRDQADataMaker::operator=(const AliTRDQADataMaker& qadm)
73 {
74   //
75   // Equal operator.
76   //
77
78   this->~AliTRDQADataMaker();
79   new(this) AliTRDQADataMaker(qadm);
80   return *this;
81
82 }
83
84 //____________________________________________________________________________ 
85 void AliTRDQADataMaker::EndOfDetectorCycle(AliQAv1::TASKINDEX task, TObjArray * list)
86 {
87   //
88   // Detector specific actions at end of cycle
89   //
90   //TStopwatch watch;
91   //watch.Start();
92   ResetEventTrigClasses();
93   //
94   //AliDebug(AliQAv1::GetQADebugLevel(), Form("EndOfCycle", "Fitting RecPoints %d", task))
95   TH1F *hist = new TH1F("fitHist", "", 200, -0.5, 199.5);
96   //
97   // RS Add a loop over species
98   for (Int_t specie = 0 ; specie < AliRecoParam::kNSpecies ; specie++) {
99     if ( !AliQAv1::Instance()->IsEventSpecieSet(specie) ) continue ; 
100     SetEventSpecie(specie);
101     //
102     for (int itc=-1;itc<GetNTrigClasses();itc++) { // RS: loop over eventual clones per trigger class
103       //
104       if (task == AliQAv1::kRECPOINTS) {
105
106         //list->Print();
107         TObjArray& arrRP = *GetRecPointsDataOfTrigClass(itc); // RS Histos matching to trigger class
108         // Rec points full chambers
109         TH2* h2tmp = (TH2*) arrRP[1];
110         if (h2tmp) {
111           for (Int_t i=0; i<540; i++) {
112             hist->Reset();
113             for(Int_t b=1; b<hist->GetXaxis()->GetNbins()-1; b++) {
114               Double_t xvalue = hist->GetBinCenter(b);
115               Int_t bin = h2tmp->FindBin(i,xvalue);
116               Double_t value =  h2tmp->GetBinContent(bin);
117               hist->SetBinContent(b, value);
118             }
119             //printf("Sum = %d %f\n", i, hist->GetSum());
120             if (hist->GetSum() < 100) continue; // chamber not present
121             //
122             hist->Fit("landau", "q0", "goff", 10, 180);
123             TF1 *fit = hist->GetFunction("landau");
124             if (arrRP[12]) ((TH1*)arrRP[12])->Fill(fit->GetParameter(1));
125             if (arrRP[13]) ((TH1*)arrRP[13])->Fill(fit->GetParameter(2));
126           }
127         }
128         //
129         // time-bin by time-bin
130         TH3* h3tmp = (TH3*) arrRP[10]; 
131         if (h3tmp) {
132           for (Int_t i=0; i<540; i++) {
133             for(Int_t j=0; j<35; j++) {
134               hist->Reset();
135               for(Int_t b=1; b<hist->GetXaxis()->GetNbins()-1; b++) {
136                 Double_t xvalue = hist->GetBinCenter(b);
137                 Int_t bin = h3tmp->FindBin(i,j,xvalue);
138                 Double_t value =  h2tmp->GetBinContent(bin);
139                 //printf("v = %f\n", value);
140                 hist->SetBinContent(b, value);
141               }
142               if (hist->GetSum() < 100) continue;
143               //printf("fitting %d %d %f\n", i, j, hist->GetSum());           
144               hist->Fit("landau", "q0", "goff", 10, 180);
145               TF1 *fit = hist->GetFunction("landau");
146               //
147               Int_t sm = i/18;
148               Int_t det = i%18;
149               TH2* h2 = (TH2*)arrRP[14+sm];
150               if (!h2) continie;
151               Int_t bin = h2->FindBin(det,j);
152               // printf("%d %d %d\n", det, j, bin);
153               h2->SetBinContent(bin, fit->GetParameter(1));
154             }
155           }
156         } // h3tmp
157       } // RESPOINTS
158     } // RS: loop over eventual clones per trigger class
159   } // loop over species
160   
161   
162   delete hist;
163   
164   // call the checker
165   AliQAChecker::Instance()->Run(AliQAv1::kTRD, task, list) ;    
166
167   //watch.Stop();
168   //watch.Print();
169 }
170
171 //____________________________________________________________________________ 
172 void AliTRDQADataMaker::InitESDs()
173 {
174   //
175   // Create ESDs histograms in ESDs subdir
176   //
177
178   const Int_t kNhist = 19;
179   TH1 *hist[kNhist];
180   Int_t histoCounter = -1 ;
181
182   hist[++histoCounter] = new TH1F("qaTRD_esd_ntracks", ":Number of tracks", 300, -0.5, 299.5);
183   hist[++histoCounter] = new TH1F("qaTRD_esd_sector", ":Sector", 18, -0.5, 17.7);
184   hist[++histoCounter] = new TH1F("qaTRD_esd_bits", ";Bits", 64, -0.5, 63.5);
185
186   const Int_t knbits = 6;
187   const char *suf[knbits] = {"TPCi", "TPCo", "TPCz", "TRDo", "TRDr", "TRDz"};
188
189   for(Int_t i=0; i<knbits; i++) {
190     hist[++histoCounter] = new TH1F(Form("qaTRD_esd_pt%s",suf[i]), ";p_{T} (GeV/c);", 50, 0, 10);
191     hist[++histoCounter] = new TH1F(Form("qaTRD_esd_trdz%s", suf[i]), ";z (cm)", 200, -400, 400); 
192   }
193
194   hist[++histoCounter] = new TH1F("qaTRD_esd_clsTRDo", "TRDo;number of clusters", 130, -0.5, 129.5);;
195   hist[++histoCounter] = new TH1F("qaTRD_esd_clsTRDr", "TRDr;number of clusters", 130, -0.5, 129.5);;
196   hist[++histoCounter] = new TH1F("qaTRD_esd_clsTRDz", "TRDz;number of clusters", 130, -0.5, 129.5);;
197   //hist[++histoCounter] = new TH1F("qaTRD_esd_clsRatio", ";cluster ratio", 100, 0., 1.3);;
198
199   hist[++histoCounter] = new TH2F("qaTRD_esd_sigMom", ";momentum (GeV/c);signal", 100, 0, 5, 200, 0, 1e3);
200
201   for(Int_t i=0; i<=histoCounter; i++) {
202     //hist[i]->Sumw2();
203     Add2ESDsList(hist[i], i);
204   }
205   //
206   ClonePerTrigClass(AliQAv1::kESDS); // this should be the last line
207 }
208
209 //____________________________________________________________________________ 
210 void AliTRDQADataMaker::InitHits()
211 {
212   //
213   // Create Hits histograms in Hits subdir
214   //
215
216   const Int_t kNhist = 4;
217   TH1F *hist[kNhist];
218
219   hist[0] = new TH1F("qaTRD_hits_det", ";Detector Id of the hit", 540, -0.5, 539.5) ; 
220
221   hist[1] = new TH1F("qaTRD_hist_Qdrift", ";Charge from tracks", 100, 0, 100);
222   hist[2] = new TH1F("qaTRD_hist_Qamp", ";Charge from TRD photon", 100, 0, 100);
223   hist[3] = new TH1F("qaTRD_hist_Qphoton", ";Charge from TRD photon", 100, 0, 100);
224
225   for(Int_t i=0; i<kNhist; i++) {
226     //hist[i]->Sumw2();
227     Add2HitsList(hist[i], i);
228   }
229   //
230   ClonePerTrigClass(AliQAv1::kHITS); // this should be the last line
231 }
232
233 //____________________________________________________________________________ 
234 void AliTRDQADataMaker::InitDigits()
235 {
236   //
237   // Create Digits histograms in Digits subdir
238   //
239
240   const Int_t kNhist = 3;
241   TH1F *hist[kNhist];
242
243   hist[0] = new TH1F("qaTRD_digits_det", ";Detector Id of the digit", 540, -0.5, 539.5);
244   hist[1] = new TH1F("qaTRD_digits_time", ";Time bin", 40, -0.5, 39.5);
245   hist[2] = new TH1F("qaTRD_digits_amp", ";Amplitude", 100, 0, 100.);
246
247   for(Int_t i=0; i<kNhist; i++) {
248     hist[i]->Sumw2();
249     Add2DigitsList(hist[i], i);
250   }
251   //
252   ClonePerTrigClass(AliQAv1::kDIGITS); // this should be the last line
253 }
254
255 //____________________________________________________________________________ 
256 void AliTRDQADataMaker::InitRecPoints()
257 {
258   //
259   // Create Reconstructed Points histograms in RecPoints subdir
260   //
261
262   const Int_t kNhist = 14 + 18;
263   TH1 *hist[kNhist];
264
265   hist[0] = new TH1F("qaTRD_recPoints_det", ";Detector ID of the cluster", 540, -0.5, 539.5);
266   hist[1] = new TH2F("qaTRD_recPoints_amp", ";Amplitude", 540, -0.5, 539, 200, -0.5, 199.5);
267   hist[2] = new TH1F("qaTRD_recPoints_npad", ";Number of Pads", 12, -0.5, 11.5);
268
269   hist[3] = new TH1F("qaTRD_recPoints_dist2", ";residuals [2pad]", 100, -1, 1);
270   hist[4] = new TH1F("qaTRD_recPoints_dist3", ";residuals [3pad]", 100, -1, 1);
271   hist[5] = new TH1F("qaTRD_recPoints_dist4", ";residuals [4pad]", 100, -1, 1);
272   hist[6] = new TH1F("qaTRD_recPoints_dist5", ";residuals [5pad]", 100, -1, 1);
273
274   hist[7] = new TH2F("qaTRD_recPoints_rowCol", ";row;col", 16, -0.5, 15.5, 145, -0.5, 144.5);
275   hist[8] = new TH1F("qaTRD_recPoints_time", ";time bin", 35, -0.5, 34.5);
276   hist[9] = new TH1F("qaTRD_recPoints_nCls", ";number of clusters", 500, -0.5, 499.5);
277
278   hist[10] = new TH3F("qaTRD_recPoints_sigTime", ";chamber;time bin;signal", 
279                       540, -0.5, 539.5, 35, -0.5, 34.5, 200, 0.5, 199.5);
280   hist[11] = new TProfile("qaTRD_recPoints_prf", ";distance;center of gravity"
281                          , 120, -0.6, 0.6, -1.2, 1.2, "");
282
283   hist[12] = new TH1F("qaTRD_recPoints_ampMPV", ";amplitude MPV", 100, 0, 100);
284   hist[13] = new TH1F("qaTRD_recPoints_ampSigma", ";amplitude Sigma", 100, 0, 100); 
285
286   for(Int_t i=0; i<18; i++) {
287     hist[14+i] = new TH2F(Form("qaTRD_recPoints_sigTime_sm%d",i), Form("sm%d;det;time bin"), 
288                         30, -0.5, 29.5, 35, -0.5, 34.5);
289     hist[14+i]->SetMinimum(20);
290     hist[14+i]->SetMaximum(40);
291   }
292
293   for(Int_t i=0; i<kNhist; i++) {
294     //hist[i]->Sumw2();
295     Add2RecPointsList(hist[i], i);
296   }
297   //
298   ClonePerTrigClass(AliQAv1::kRECPOINTS); // this should be the last line
299 }
300
301 //____________________________________________________________________________ 
302 void AliTRDQADataMaker::InitRaws()
303 {
304   //
305   // create Raws histograms in Raws subdir
306   //
307
308   const Int_t kSM = 18;
309   //const Int_t kNCh = 540;
310   const Int_t kNhist = 4+kSM;
311   TH1F *hist[kNhist];
312
313   // four histograms to be published
314   hist[0] = new TH1F("qaTRD_raws_det", ";detector", 540, -0.5, 539.5);
315   hist[1] = new TH1F("qaTRD_raws_sig", ";signal", 100, -0.5, 99.5);
316   hist[2] = new TH1F("qaTRD_raws_timeBin", ";time bin", 40, -0.5, 39.5); 
317   hist[3] = new TH1F("qaTRD_raws_smId", ";supermodule", 18, -0.5, 17.5);
318   //
319   
320   // one double per MCM (not published)
321   const Int_t kNMCM = 30 * 8 * 16;
322   for(Int_t i=0; i<kSM; i++)
323     hist[4+i] = new TH1F(Form("qaTRD_raws_sm%d",i),"",kNMCM, -0.5, kNMCM-0.5); 
324   
325   // register
326   for(Int_t i=0; i<kNhist; i++) {
327     //hist[i]->Sumw2();
328     Add2RawsList(hist[i], i);
329   }
330   //
331   ClonePerTrigClass(AliQAv1::kRAWS); // this should be the last line
332 }
333
334 //____________________________________________________________________________ 
335 void AliTRDQADataMaker::InitSDigits()
336 {
337   //
338   // Create SDigits histograms in SDigits subdir
339   //
340
341   const Int_t kNhist = 3;
342   TH1F *hist[kNhist];
343
344   hist[0] = new TH1F("qaTRD_sdigits_det", ";Detector Id of the digit", 540, -0.5, 539.5);
345   hist[1] = new TH1F("qaTRD_sdigits_time", ";Time bin", 40, -0.5, 39.5);
346   hist[2] = new TH1F("qaTRD_sdigits_amp", ";Amplitude", 100, 0, 1e7);
347
348   for(Int_t i=0; i<kNhist; i++) {
349     hist[i]->Sumw2();
350     Add2SDigitsList(hist[i], i);
351   }
352   //
353   ClonePerTrigClass(AliQAv1::kSDIGITS); // this should be the last line
354 }
355
356 //____________________________________________________________________________
357 void AliTRDQADataMaker::MakeESDs(AliESDEvent * const esd)
358 {
359   //
360   // Make QA data from ESDs
361   //
362
363   Int_t nTracks = esd->GetNumberOfTracks();
364   FillESDsData(0,nTracks);
365
366   // track loop
367   for (Int_t i=0; i<nTracks; i++) {
368
369     AliESDtrack *track = esd->GetTrack(i);
370     const AliExternalTrackParam *paramOut = track->GetOuterParam();
371     const AliExternalTrackParam *paramIn = track->GetInnerParam();
372
373     // long track ..
374     if (!paramIn) continue;
375     if (!paramOut) continue;
376
377     // not a kink
378     if (track->GetKinkIndex(0) > 0) continue; 
379
380     Double_t extZ = GetExtZ(paramIn);
381     if (TMath::Abs(extZ) > 320) continue; // acceptance cut
382
383     // .. in the acceptance
384     Int_t sector = GetSector(paramOut->GetAlpha());
385
386     UInt_t u = 1;
387     UInt_t status = track->GetStatus();
388     for(Int_t bit=0; bit<32; bit++) 
389       if (u<<bit & status) FillESDsData(2,bit);
390
391     const Int_t knbits = 6; 
392     Int_t bit[6] = {0,0,0,0,0,0};    
393     bit[0] = status & AliESDtrack::kTPCin;
394     bit[1] = status & AliESDtrack::kTPCout;
395     bit[2] = (status & AliESDtrack::kTPCout) && !(status & AliESDtrack::kTRDout);
396     bit[3] = status & AliESDtrack::kTRDout;
397     bit[4] = status & AliESDtrack::kTRDrefit;
398     bit[5] = (status & AliESDtrack::kTRDout) && !(status & AliESDtrack::kTRDrefit);
399
400     // transverse momentum
401     //const Double_t *val = paramOut->GetParameter(); // parameters at the Outer plane
402     Double_t pt = paramOut->Pt(); //1./TMath::Abs(val[4]);
403
404     for(Int_t b=0; b<knbits; b++) {
405       if (bit[b]) {
406         FillESDsData(2*b+3,pt); 
407         FillESDsData(2*b+4,extZ);
408       }
409     }
410
411     // clusters
412     for(Int_t b=0; b<3; b++) 
413       if (bit[3+b]) FillESDsData(b+15,track->GetTRDncls());
414
415     // refitted only
416     if (!bit[4]) continue;
417
418     //fQuality->Fill(track->GetTRDQuality());
419     //fBudget->Fill(track->GetTRDBudget());
420     //fSignal->Fill(track->GetTRDsignal());
421         
422     FillESDsData(18,track->GetP(), track->GetTRDsignal());
423     FillESDsData(1,sector);
424
425     /*
426     // PID only
427     if (status & AliESDtrack::kTRDpid) {
428
429       for(Int_t l=0; l<6; l++) fTime->Fill(track->GetTRDTimBin(l));
430
431       // fill pid histograms
432       Double_t trdr0 = 0; //, tpcr0 = 0;
433       Int_t trdBestPid = 5; //, tpcBestPid = 5;  // charged
434       const Double_t kminPidValue = 0.9;
435
436       //Double_t pp[5];
437       //track->GetTPCpid(pp); // ESD inconsequence
438
439       for(Int_t pid=0; pid<5; pid++) {
440         
441         trdr0 += track->GetTRDpid(pid);
442         //tpcr0 += pp[pid];
443         
444         fTrdPID[pid]->Fill(track->GetTRDpid(pid));
445         //fTpcPID[pid]->Fill(pp[pid]);
446         
447         if (track->GetTRDpid(pid) > kminPidValue) trdBestPid = pid;
448         //if (pp[pid] > kminPidValue) tpcBestPid = pid;
449       }
450
451       fTrdPID[5]->Fill(trdr0); // check unitarity
452       fTrdSigMomPID[trdBestPid]->Fill(track->GetP(), track->GetTRDsignal());
453
454       //fTpcPID[5]->Fill(tpcr0); // check unitarity
455       //fTpcSigMomPID[tpcBestPid]->Fill(track->GetP(), track->GetTPCsignal());
456     }
457     */
458
459   }
460   //
461   IncEvCountCycleESDs();
462   IncEvCountTotalESDs();
463   //
464 }
465
466 //______________________________________________________________________________
467 Int_t AliTRDQADataMaker::GetSector(const Double_t alpha) const 
468 {
469   //
470   // Gets the sector number 
471   //
472
473   Double_t size = TMath::DegToRad() * 20.; // shall use TRDgeo
474   Int_t sector = (Int_t)((alpha + TMath::Pi())/size);
475   return sector;
476
477 }
478
479 //______________________________________________________________________________
480 Double_t AliTRDQADataMaker::GetExtZ(const AliExternalTrackParam *in) const 
481 {
482   //
483   // Returns the Z position at the entry to TRD
484   // using parameters from the TPC in
485   //
486
487   const Double_t kX0 = 300;
488
489   Double_t x = in->GetX();
490   const Double_t *par = in->GetParameter();
491   Double_t theta = par[3];
492   Double_t z = in->GetZ();
493
494   Double_t zz = z + (kX0-x) * TMath::Tan(theta);
495   return zz;
496
497 }
498
499 //____________________________________________________________________________
500 void AliTRDQADataMaker::MakeHits(TClonesArray * const hits)
501 {
502   //
503   // Make QA data from Hits
504   //
505
506   TIter next(hits); 
507   AliTRDhit * hit; 
508
509   while ( (hit = dynamic_cast<AliTRDhit *>(next())) ) {
510     FillHitsData(0,hit->GetDetector());
511     Double_t q = TMath::Abs(hit->GetCharge());
512
513     if (hit->FromDrift()) FillHitsData(1,q);
514     if (hit->FromAmplification()) FillHitsData(2,q);
515     if (hit->FromTRphoton()) FillHitsData(3,q);
516   }
517
518 }
519
520 //____________________________________________________________________________
521 void AliTRDQADataMaker::MakeHits(TTree * hitTree)
522 {
523   //
524   // Make QA data from Hits
525   //
526
527   if (!CheckPointer(hitTree, "TRD hits tree")) return;
528
529   TBranch *branch = hitTree->GetBranch("TRD");
530   if (!CheckPointer(branch, "TRD hits branch")) return;
531
532   Int_t nhits = (Int_t)(hitTree->GetTotBytes()/sizeof(AliTRDhit));
533   TClonesArray *hits = new TClonesArray("AliTRDhit", nhits+1000);
534   TClonesArray *tmp = new TClonesArray("AliTRDhit", 1000);
535   branch->SetAddress(&tmp);
536
537   Int_t index = 0;
538   Int_t nEntries = (Int_t)branch->GetEntries();
539   for(Int_t i = 0; i < nEntries; i++) {
540     branch->GetEntry(i);
541     Int_t nHits = (Int_t)tmp->GetEntries();
542     for(Int_t j=0; j<nHits; j++) {
543       AliTRDhit *hit = (AliTRDhit*)tmp->At(j);
544       new((*hits)[index++]) AliTRDhit(*hit);
545     }
546   }
547
548   tmp->Delete();
549   delete tmp;
550   MakeHits(hits);
551   //
552   IncEvCountCycleHits();
553   IncEvCountTotalHits();
554   //
555 }
556
557 //____________________________________________________________________________
558 void AliTRDQADataMaker::MakeDigits(TClonesArray * const digits)
559 {
560   //
561   // Makes data from Digits
562   //
563
564   TIter next(digits) ; 
565   AliTRDdigit * digit ; 
566   while ( (digit = dynamic_cast<AliTRDdigit *>(next())) ) {
567     FillDigitsData(0,digit->GetDetector());
568     FillDigitsData(1,digit->GetTime());
569     FillDigitsData(2,digit->GetAmp());
570   }
571
572 }
573
574 //____________________________________________________________________________
575 void AliTRDQADataMaker::MakeDigits(TTree * digits)
576 {
577   //
578   // Makes data from digits tree
579   //
580
581   AliTRDdigitsManager *digitsManager = new AliTRDdigitsManager();
582   digitsManager->CreateArrays();
583   digitsManager->ReadDigits(digits);
584
585   for (Int_t i = 0; i < AliTRDgeometry::kNdet; i++) {
586
587     AliTRDarrayADC *digitsIn = (AliTRDarrayADC *) digitsManager->GetDigits(i);      
588
589     // This is to take care of switched off super modules
590     if (digitsIn->GetNtime() == 0) continue;
591
592     digitsIn->Expand(); 
593
594     //AliTRDSignalIndex* indexes = digitsManager->GetIndexes(i);
595     //if (indexes->IsAllocated() == kFALSE) digitsManager->BuildIndexes(i);
596
597     Int_t nRows = digitsIn->GetNrow();
598     Int_t nCols = digitsIn->GetNcol();
599     Int_t nTbins = digitsIn->GetNtime();
600
601     for(Int_t row = 0; row < nRows; row++) 
602       for(Int_t col = 0; col < nCols; col++) 
603         for(Int_t time = 0; time < nTbins; time++) 
604           {
605             Float_t signal = digitsIn->GetData(row,col,time);
606             FillDigitsData(0,i);
607             FillDigitsData(1,time);
608             FillDigitsData(2,signal);
609           }
610
611     //delete digitsIn;
612   }
613
614   delete digitsManager;
615   //
616   IncEvCountCycleDigits();
617   IncEvCountTotalDigits();
618   //
619 }
620
621 //____________________________________________________________________________
622 void AliTRDQADataMaker::MakeSDigits(TClonesArray * const sdigits)
623 {
624   //
625   // Makes data from Digits
626   //
627
628   TIter next(sdigits) ; 
629   AliTRDdigit * digit ; 
630   while ( (digit = dynamic_cast<AliTRDdigit *>(next())) ) {
631     FillDigitsData(0,digit->GetDetector());
632     FillDigitsData(1,digit->GetTime());
633     FillDigitsData(2,digit->GetAmp());
634   }
635
636 }
637
638 //____________________________________________________________________________
639 void AliTRDQADataMaker::MakeSDigits(TTree * digits)
640 {
641   //
642   // Makes data from SDigits
643   //
644
645   AliTRDdigitsManager *digitsManager = new AliTRDdigitsManager();
646   digitsManager->CreateArrays();
647   digitsManager->ReadDigits(digits);
648
649   for (Int_t i = 0; i < AliTRDgeometry::kNdet; i++) 
650     {
651     AliTRDarrayADC *digitsIn = (AliTRDarrayADC *) digitsManager->GetDigits(i);      
652
653     // This is to take care of switched off super modules
654     if (digitsIn->GetNtime() == 0) continue;
655
656     digitsIn->Expand();  
657
658     //AliTRDSignalIndex* indexes = digitsManager->GetIndexes(i);
659     //if (indexes->IsAllocated() == kFALSE) digitsManager->BuildIndexes(i);
660
661     Int_t nRows = digitsIn->GetNrow();
662     Int_t nCols = digitsIn->GetNcol();
663     Int_t nTbins = digitsIn->GetNtime();
664
665     for(Int_t row = 0; row < nRows; row++) 
666       for(Int_t col = 0; col < nCols; col++) 
667         for(Int_t time = 0; time < nTbins; time++) 
668           {
669             Float_t signal = digitsIn->GetData(row,col,time);
670             if (signal <= 0) continue;
671             FillSDigitsData(0,i);
672             FillSDigitsData(1,time);
673             FillSDigitsData(2,signal);
674           }
675     
676     // delete digitsIn;
677   }
678
679   delete digitsManager;
680   //
681   IncEvCountCycleSDigits();
682   IncEvCountTotalSDigits();
683   //
684 }
685
686 //____________________________________________________________________________
687 void AliTRDQADataMaker::MakeRaws(AliRawReader * const rawReader)
688 {
689   //
690   // Makes QA data from raw data
691   //
692
693   // 157
694   // T9 -- T10
695
696   //const Int_t kSM = 18;
697   //const Int_t kROC = 30;
698   const Int_t kROB = 8;
699   //const Int_t kLayer = 6;
700   //const Int_t kStack = 5;
701   const Int_t kMCM = 16;
702   //  const Int_t kADC = 22;
703
704   AliTRDrawStream raw(rawReader);
705   raw.SetRawVersion(3);
706   raw.Init();
707
708   while (raw.Next()) {
709
710     FillRawsData(0,raw.GetDet());
711
712     // possibly needs changes with the new reader !!
713     Int_t *sig = raw.GetSignals();
714     for(Int_t i=0; i<3; i++) FillRawsData(1,sig[i]);
715     // ---
716
717     FillRawsData(2,raw.GetTimeBin());
718
719     // calculate the index;
720     Int_t sm = raw.GetSM();
721     Int_t roc = raw.GetROC();
722     Int_t rob = raw.GetROB();
723     Int_t mcm = raw.GetMCM();
724     //Int_t adc = raw.GetADC();
725
726     //Int_t index = roc * (kROB*kMCM*kADC) + rob * (kMCM*kADC) + mcm * kADC + adc;
727     Int_t  index = roc * (kROB*kMCM) + rob * kMCM + mcm;
728     FillRawsData(3,sm);
729     FillRawsData(4+sm,index);
730   }
731   //
732   IncEvCountCycleRaws();
733   IncEvCountTotalRaws();
734   //
735 }
736
737 //____________________________________________________________________________
738 void AliTRDQADataMaker::MakeRecPoints(TTree * clustersTree)
739 {
740   //  
741   // Makes data from RecPoints
742   // 
743
744   Int_t nsize = Int_t(clustersTree->GetTotBytes() / (sizeof(AliTRDcluster))); 
745   TObjArray *clusterArray = new TObjArray(nsize+1000); 
746
747   TBranch *branch = clustersTree->GetBranch("TRDcluster");
748   if (!branch) {
749     AliError("Can't get the branch !");
750     return;
751   }
752   branch->SetAddress(&clusterArray); 
753
754   // Loop through all entries in the tree
755   Int_t nEntries   = (Int_t) clustersTree->GetEntries();
756   Int_t nbytes     = 0;
757   AliTRDcluster *c = 0;
758   Int_t nDet[540];
759   for (Int_t i=0; i<540; i++) nDet[i] = 0;
760
761   for (Int_t iEntry = 0; iEntry < nEntries; iEntry++) {    
762
763     // Import the tree
764     nbytes += clustersTree->GetEvent(iEntry);  
765
766     // Get the number of points in the detector
767     Int_t nCluster = clusterArray->GetEntriesFast();  
768
769     // Loop through all TRD digits
770     for (Int_t iCluster = 0; iCluster < nCluster; iCluster++) { 
771       c = (AliTRDcluster *) clusterArray->UncheckedAt(iCluster);
772
773       Int_t iDet = c->GetDetector();
774       nDet[iDet]++;
775       FillRecPointsData(0,iDet);
776       FillRecPointsData(1,iDet, c->GetQ());
777       FillRecPointsData(2,c->GetNPads());
778       if (c->GetNPads() < 6)
779         FillRecPointsData(1+c->GetNPads(),c->GetCenter());
780
781       //if (c->GetPadTime() < 5)
782       FillRecPointsData(7,c->GetPadRow(), c->GetPadCol());
783       FillRecPointsData(8,c->GetPadTime());
784
785       TObjArray *hists = GetMatchingRecPointsData(10); //RS no alias for 3d histo filling, to directly
786       for (int ih=hists->GetEntriesFast();ih--;) ((TH3F*)hists->UncheckedAt(ih))->Fill(iDet, c->GetPadTime(), c->GetQ());
787
788       // PRF for 2pad
789       //if (c->GetNPads() == 2) {
790       Short_t *sig = c->GetSignals();
791       Double_t frac = -10;
792
793       if (sig[0] == 0 && sig[1] == 0 && sig[2] == 0 && sig[5] == 0 && sig[6] == 0) 
794         frac = 1. * sig[4] / (sig[3] + sig[4]);
795
796       if (sig[0] == 0 && sig[1] == 0 && sig[4] == 0 && sig[5] == 0 && sig[6] == 0)
797         frac = -1. * sig[2] / (sig[2] + sig[3]);
798
799       if (frac > -10)  FillRecPointsData(11,c->GetCenter(), frac);
800         
801       //}
802     }
803   }
804
805   for(Int_t i=0; i<540; i++) 
806     if (nDet[i] > 0) FillRecPointsData(9,nDet[i]);
807
808   delete clusterArray;
809
810 }
811
812 //____________________________________________________________________________ 
813 void AliTRDQADataMaker::StartOfDetectorCycle()
814 {
815   //
816   // Detector specific actions at start of cycle
817   //
818
819 }
820
821 //__________________________________________________________________________
822 Int_t AliTRDQADataMaker::CheckPointer(TObject * const obj, const char *name) 
823 {
824   //
825   // Checks initialization of pointers
826   //
827
828   if (!obj) AliWarning(Form("null pointer: %s", name));
829   return !!obj;
830
831 }