1 /**************************************************************************
2 * Copyright(c) 1998-2007, ALICE Experiment at CERN, All rights reserved. *
4 * Author: The ALICE Off-line Project. *
5 * Contributors are mentioned in the code where appropriate. *
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 **************************************************************************/
20 Based on AliPHOSQADataMaker
21 Produces the data needed to calculate the quality assurance.
22 All data must be mergeable objects.
23 P. Christiansen, Lund, January 2008
29 We have chosen to have the histograms as non-persistent meber to
30 allow better debugging. In the copy constructor we then have to
31 assign the pointers to the existing histograms in the copied
32 list. This have been implemented but not tested.
35 #include "AliTPCQADataMakerSim.h"
37 // --- ROOT system ---
39 // --- Standard library ---
41 // --- AliRoot header files ---
42 #include "AliQAChecker.h"
45 #include "AliSimDigits.h"
48 ClassImp(AliTPCQADataMakerSim)
50 //____________________________________________________________________________
51 AliTPCQADataMakerSim::AliTPCQADataMakerSim() :
52 AliQADataMakerSim(AliQAv1::GetDetName(AliQAv1::kTPC),
53 "TPC Sim Quality Assurance Data Maker")
58 //____________________________________________________________________________
59 AliTPCQADataMakerSim::AliTPCQADataMakerSim(const AliTPCQADataMakerSim& qadm) :
63 SetName((const char*)qadm.GetName()) ;
64 SetTitle((const char*)qadm.GetTitle());
67 // Associate class histogram objects to the copies in the list
68 // Could also be done with the indexes
72 //__________________________________________________________________
73 AliTPCQADataMakerSim& AliTPCQADataMakerSim::operator = (const AliTPCQADataMakerSim& qadm )
76 this->~AliTPCQADataMakerSim();
77 new(this) AliTPCQADataMakerSim(qadm);
81 //____________________________________________________________________________
82 void AliTPCQADataMakerSim::EndOfDetectorCycle(AliQAv1::TASKINDEX_t task, TObjArray ** list)
84 //Detector specific actions at end of cycle
86 AliQAChecker::Instance()->Run(AliQAv1::kTPC, task, list) ;
89 //____________________________________________________________________________
90 void AliTPCQADataMakerSim::InitDigits()
92 const Bool_t expert = kTRUE ;
93 const Bool_t image = kTRUE ;
94 TH1F * histDigitsADC =
95 new TH1F("hDigitsADC", "Digit ADC distribution; ADC; Counts",
97 histDigitsADC->Sumw2();
98 Add2DigitsList(histDigitsADC, kDigitsADC, !expert, image);
100 ClonePerTrigClass(AliQAv1::kDIGITS); // this should be the last line
103 //____________________________________________________________________________
104 void AliTPCQADataMakerSim::InitHits()
106 const Bool_t expert = kTRUE ;
107 const Bool_t image = kTRUE ;
108 TH1F * histHitsNhits =
109 new TH1F("hHitsNhits", "Interactions per track in the TPC volume; Number of interactions; Counts",
111 histHitsNhits->Sumw2();
112 Add2HitsList(histHitsNhits, kNhits, !expert, image);
114 TH1F * histHitsElectrons =
115 new TH1F("hHitsElectrons", "Electrons per interaction; Electrons; Counts",
117 histHitsElectrons->Sumw2();
118 Add2HitsList(histHitsElectrons, kElectrons, !expert, image);
120 TH1F * histHitsRadius =
121 new TH1F("hHitsRadius", "Position of interaction; Radius; Counts",
123 histHitsRadius->Sumw2();
124 Add2HitsList(histHitsRadius, kRadius, !expert, image);
126 TH1F * histHitsPrimPerCm =
127 new TH1F("hHitsPrimPerCm", "Primaries per cm; Primaries; Counts",
129 histHitsPrimPerCm->Sumw2();
130 Add2HitsList(histHitsPrimPerCm, kPrimPerCm, !expert, image);
132 TH1F * histHitsElectronsPerCm =
133 new TH1F("hHitsElectronsPerCm", "Electrons per cm; Electrons; Counts",
135 histHitsElectronsPerCm->Sumw2();
136 Add2HitsList(histHitsElectronsPerCm, kElectronsPerCm, !expert, image);
138 ClonePerTrigClass(AliQAv1::kHITS); // this should be the last line
141 //____________________________________________________________________________
142 void AliTPCQADataMakerSim::MakeDigits(TTree* digitTree)
145 TBranch* branch = digitTree->GetBranch("Segment");
146 AliSimDigits* digArray = 0;
147 branch->SetAddress(&digArray);
149 Int_t nEntries = Int_t(digitTree->GetEntries());
151 for (Int_t n = 0; n < nEntries; n++) {
153 digitTree->GetEvent(n);
155 if (digArray->First())
157 Float_t dig = digArray->CurrentDigit();
159 FillDigitsData(kDigitsADC,dig);
160 } while (digArray->Next());
163 IncEvCountCycleDigits();
164 IncEvCountTotalDigits();
168 //____________________________________________________________________________
169 void AliTPCQADataMakerSim::MakeHits(TTree * hitTree)
171 // make QA data from Hit Tree
173 const Int_t nTracks = hitTree->GetEntries();
174 TBranch* branch = hitTree->GetBranch("TPC2");
175 AliTPCv2* tpc = (AliTPCv2*)gAlice->GetDetector("TPC");
180 for(Int_t n = 0; n < nTracks; n++){
184 AliTPChit* tpcHit = (AliTPChit*)tpc->FirstHit(-1);
189 Float_t xold = tpcHit->X();
190 Float_t yold = tpcHit->Y();
191 Float_t zold = tpcHit->Z();
192 Float_t radiusOld = TMath::Sqrt(xold*xold + yold*yold);
193 Int_t trackOld = tpcHit->GetTrack();
198 Float_t x = tpcHit->X();
199 Float_t y = tpcHit->Y();
200 Float_t z = tpcHit->Z();
201 Float_t radius = TMath::Sqrt(x*x + y*y);
203 if(radius>50) { // Skip hits at interaction point
207 Int_t trackNo = tpcHit->GetTrack();
209 FillHitsData(kElectrons,tpcHit->fQ);
210 FillHitsData(kRadius,radius);
212 if(trackNo==trackOld) { // if the same track
214 // find the new distance
215 dist += TMath::Sqrt((x-xold)*(x-xold) + (y-yold)*(y-yold) +
217 if(dist<1.){ // add data to this 1 cm step
221 } else{ // Fill the histograms normalized to per cm
224 // cout << radius << ", " << radiusOld << ", " << dist << endl;
226 FillHitsData(kPrimPerCm,(Float_t)nprim);
227 FillHitsData(kElectronsPerCm,q);
233 } else { // reset for new track
245 trackOld = tpcHit->GetTrack();
247 tpcHit = (AliTPChit*) tpc->NextHit();
251 FillHitsData(kNhits,nHits);
254 IncEvCountCycleHits();
255 IncEvCountTotalHits();