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"
46 #include "AliSimDigits.h"
48 ClassImp(AliTPCQADataMakerSim)
50 //____________________________________________________________________________
51 AliTPCQADataMakerSim::AliTPCQADataMakerSim() :
52 AliQADataMakerSim(AliQA::GetDetName(AliQA::kTPC),
53 "TPC Sim Quality Assurance Data Maker"),
55 fHistHitsNhits(0), fHistHitsElectrons(0), fHistHitsRadius(0),
56 fHistHitsPrimPerCm(0), fHistHitsElectronsPerCm(0)
61 //____________________________________________________________________________
62 AliTPCQADataMakerSim::AliTPCQADataMakerSim(const AliTPCQADataMakerSim& qadm) :
66 SetName((const char*)qadm.GetName()) ;
67 SetTitle((const char*)qadm.GetTitle());
70 // Associate class histogram objects to the copies in the list
71 // Could also be done with the indexes
73 fHistDigitsADC = (TH1F*)fDigitsQAList->FindObject("hDigitsADC");
75 fHistHitsNhits = (TH1F*)fHitsQAList->FindObject("hHitsNhits");
76 fHistHitsElectrons = (TH1F*)fHitsQAList->FindObject("hHitsElectrons");
77 fHistHitsRadius = (TH1F*)fHitsQAList->FindObject("hHitsRadius");
78 fHistHitsPrimPerCm = (TH1F*)fHitsQAList->FindObject("hHitsPrimPerCm");
79 fHistHitsElectronsPerCm = (TH1F*)fHitsQAList->FindObject("hHitsElectronsPerCm");
82 //__________________________________________________________________
83 AliTPCQADataMakerSim& AliTPCQADataMakerSim::operator = (const AliTPCQADataMakerSim& qadm )
86 this->~AliTPCQADataMakerSim();
87 new(this) AliTPCQADataMakerSim(qadm);
91 //____________________________________________________________________________
92 void AliTPCQADataMakerSim::EndOfDetectorCycle(AliQA::TASKINDEX task, TObjArray * list)
94 //Detector specific actions at end of cycle
96 AliQAChecker::Instance()->Run(AliQA::kTPC, task, list) ;
99 //____________________________________________________________________________
100 void AliTPCQADataMakerSim::InitDigits()
103 new TH1F("hDigitsADC", "Digit ADC distribution; ADC; Counts",
105 fHistDigitsADC->Sumw2();
106 Add2DigitsList(fHistDigitsADC, 0);
109 //____________________________________________________________________________
110 void AliTPCQADataMakerSim::InitHits()
113 new TH1F("hHitsNhits", "Interactions per primary track in the TPC volume; Number of primary interactions; Counts",
115 fHistHitsNhits->Sumw2();
116 Add2HitsList(fHistHitsNhits, 0);
119 new TH1F("hHitsElectrons", "Electrons per interaction (primaries); Electrons; Counts",
121 fHistHitsElectrons->Sumw2();
122 Add2HitsList(fHistHitsElectrons, 1);
125 new TH1F("hHitsRadius", "Position of interaction (Primary tracks only); Radius; Counts",
127 fHistHitsRadius->Sumw2();
128 Add2HitsList(fHistHitsRadius, 2);
131 new TH1F("hHitsPrimPerCm", "Primaries per cm (Primary tracks only); Primaries; Counts",
133 fHistHitsPrimPerCm->Sumw2();
134 Add2HitsList(fHistHitsPrimPerCm, 3);
136 fHistHitsElectronsPerCm =
137 new TH1F("hHitsElectronsPerCm", "Electrons per cm (Primary tracks only); Electrons; Counts",
139 fHistHitsElectronsPerCm->Sumw2();
140 Add2HitsList(fHistHitsElectronsPerCm, 4);
143 //____________________________________________________________________________
144 void AliTPCQADataMakerSim::MakeDigits(TTree* digitTree)
146 TBranch* branch = digitTree->GetBranch("Segment");
147 AliSimDigits* digArray = 0;
148 branch->SetAddress(&digArray);
150 Int_t nEntries = Int_t(digitTree->GetEntries());
152 for (Int_t n = 0; n < nEntries; n++) {
154 digitTree->GetEvent(n);
156 if (digArray->First())
159 Float_t dig = digArray->CurrentDigit();
161 fHistDigitsADC->Fill(dig);
162 } while (digArray->Next());
166 //____________________________________________________________________________
167 void AliTPCQADataMakerSim::MakeHits(TTree * hitTree)
169 // make QA data from Hit Tree
170 const Int_t nTracks = hitTree->GetEntries();
171 TBranch* branch = hitTree->GetBranch("TPC2");
172 AliTPCv2* tpc = (AliTPCv2*)gAlice->GetDetector("TPC");
177 for(Int_t n = 0; n < nTracks; n++){
181 AliTPChit* tpcHit = (AliTPChit*)tpc->FirstHit(-1);
186 Float_t xold = tpcHit->X();
187 Float_t yold = tpcHit->Y();
188 Float_t zold = tpcHit->Z();
189 Float_t radiusOld = TMath::Sqrt(xold*xold + yold*yold);
192 Float_t x = tpcHit->X();
193 Float_t y = tpcHit->Y();
194 Float_t z = tpcHit->Z();
195 Float_t radius = TMath::Sqrt(x*x + y*y);
197 if(radius>50) { // Skip hits at interaction point
201 Int_t trackNo = tpcHit->GetTrack();
203 if(trackNo==n) { // primary track
205 fHistHitsElectrons->Fill(tpcHit->fQ);
206 fHistHitsRadius->Fill(radius);
208 // find the new distance
209 dist += TMath::Sqrt((x-xold)*(x-xold) + (y-yold)*(y-yold) +
212 // add data to this 1 cm step
217 // Fill the histograms normalized to per cm
220 cout << radius << ", " << radiusOld << ", " << dist << endl;
222 fHistHitsPrimPerCm->Fill((Float_t)nprim);
223 fHistHitsElectronsPerCm->Fill(q);
237 tpcHit = (AliTPChit*) tpc->NextHit();
240 fHistHitsNhits->Fill(nHits);
244 //____________________________________________________________________________
245 void AliTPCQADataMakerSim::MakeSDigits(TTree* sdigitTree)
247 // This method is currently not supproted for the TPC as the
248 // summable digits are never created
249 AliInfo("AliTPCQADataMakerSim::MakeSDigits() method not supprted");