1 /**************************************************************************
\r
2 * Copyright(c) 2007-2009, ALICE Experiment at CERN, All rights reserved. *
\r
4 * Author: The ALICE Off-line Project. *
\r
5 * Contributors are mentioned in the code where appropriate. *
\r
7 * Permission to use, copy, modify and distribute this software and its *
\r
8 * documentation strictly for non-commercial purposes is hereby granted *
\r
9 * without fee, provided that the above copyright notice appears in all *
\r
10 * copies and that both the copyright notice and this permission notice *
\r
11 * appear in the supporting documentation. The authors make no claims *
\r
12 * about the suitability of this software for any purpose. It is *
\r
13 * provided "as is" without express or implied warranty. *
\r
14 **************************************************************************/
\r
18 // *************************************************************
\r
19 // Checks the quality assurance
\r
20 // by comparing with reference data
\r
21 // contained in a DB
\r
22 // -------------------------------------------------------------
\r
23 // W. Ferrarese + P. Cerello Feb 2008
\r
26 // --- ROOT system ---
\r
27 #include <TProfile2D.h>
\r
29 #include <TBranch.h>
\r
33 // --- Standard library ---
\r
35 // --- AliRoot header files ---
\r
36 #include "AliITSQASDDDataMakerRec.h"
\r
39 #include "AliQAChecker.h"
\r
40 #include "AliRawReader.h"
\r
41 #include "AliITSRawStreamSDD.h"
\r
42 #include "AliITSRecPoint.h"
\r
43 #include "AliITSgeomTGeo.h"
\r
45 #include "AliCDBManager.h"
\r
46 #include "AliCDBStorage.h"
\r
47 #include "AliCDBEntry.h"
\r
50 ClassImp(AliITSQASDDDataMakerRec)
\r
52 //____________________________________________________________________________
\r
53 AliITSQASDDDataMakerRec::AliITSQASDDDataMakerRec(AliITSQADataMakerRec *aliITSQADataMakerRec, Bool_t kMode, Short_t ldc) :
\r
55 fAliITSQADataMakerRec(aliITSQADataMakerRec),
\r
63 //ctor used to discriminate OnLine-Offline analysis
\r
64 if(fLDC < 0 || fLDC > 4) {
\r
65 AliError("Error: LDC number out of range; return\n");
\r
69 //____________________________________________________________________________
\r
70 AliITSQASDDDataMakerRec::AliITSQASDDDataMakerRec(const AliITSQASDDDataMakerRec& qadm) :
\r
72 fAliITSQADataMakerRec(qadm.fAliITSQADataMakerRec),
\r
73 fkOnline(qadm.fkOnline),
\r
75 fSDDhTask(qadm.fSDDhTask),
\r
76 fGenOffset(qadm.fGenOffset),
\r
81 fAliITSQADataMakerRec->SetName((const char*)qadm.fAliITSQADataMakerRec->GetName()) ;
\r
82 fAliITSQADataMakerRec->SetTitle((const char*)qadm.fAliITSQADataMakerRec->GetTitle());
\r
85 //____________________________________________________________________________
\r
86 AliITSQASDDDataMakerRec::~AliITSQASDDDataMakerRec(){
\r
88 // if(fDDLModuleMap) delete fDDLModuleMap;
\r
90 //__________________________________________________________________
\r
91 AliITSQASDDDataMakerRec& AliITSQASDDDataMakerRec::operator = (const AliITSQASDDDataMakerRec& qac )
\r
94 this->~AliITSQASDDDataMakerRec();
\r
95 new(this) AliITSQASDDDataMakerRec(qac);
\r
99 //____________________________________________________________________________
\r
100 void AliITSQASDDDataMakerRec::StartOfDetectorCycle()
\r
102 //Detector specific actions at start of cycle
\r
103 AliDebug(1,"AliITSQADM::Start of SDD Cycle\n");
\r
106 //____________________________________________________________________________
\r
107 void AliITSQASDDDataMakerRec::EndOfDetectorCycle(AliQA::TASKINDEX_t /*task*/, TObjArray* /*list*/)
\r
109 // launch the QA checking
\r
110 AliDebug(1,"AliITSDM instantiates checker with Run(AliQA::kITS, task, list)\n");
\r
113 //____________________________________________________________________________
\r
114 void AliITSQASDDDataMakerRec::InitRaws()
\r
116 // Initialization for RAW data - SDD -
\r
117 fGenOffset = (fAliITSQADataMakerRec->fRawsQAList)->GetEntries();
\r
119 AliCDBEntry *ddlMapSDD = AliCDBManager::Instance()->Get("ITS/Calib/DDLMapSDD");
\r
120 Bool_t cacheStatus = AliCDBManager::Instance()->GetCacheFlag();
\r
123 AliError("Calibration object retrieval failed! SDD will not be processed");
\r
124 fDDLModuleMap = NULL;
\r
127 fDDLModuleMap = (AliITSDDLModuleMapSDD*)ddlMapSDD->GetObject();
\r
128 if(!cacheStatus)ddlMapSDD->SetObject(NULL);
\r
129 ddlMapSDD->SetOwner(kTRUE);
\r
130 if(!cacheStatus)delete ddlMapSDD;
\r
132 Int_t lay, lad, det;
\r
133 Int_t LAY = -1; //, LAD = -1;
\r
135 Int_t indexlast = 0;
\r
138 if(fLDC == 1 || fLDC == 2) LAY = 2;
\r
139 if(fLDC == 3 || fLDC == 4) LAY = 3;
\r
142 AliInfo("Book Online Histograms for SDD\n");
\r
145 AliInfo("Book Offline Histograms for SDD\n ");
\r
147 TH1D *h0 = new TH1D("ModPattern","HW Modules pattern",fgknSDDmodules,-0.5,259.5);
\r
148 h0->GetXaxis()->SetTitle("Module Number");
\r
149 h0->GetYaxis()->SetTitle("Counts");
\r
150 fAliITSQADataMakerRec->Add2RawsList((new TH1D(*h0)),0+fGenOffset);
\r
153 if(fLDC==0 || fLDC==1 || fLDC==2){
\r
154 TH1D *h1 = new TH1D("LadPatternL3","Ladder pattern L3",14,0.5,14.5);
\r
155 h1->GetXaxis()->SetTitle("Ladder Number on Lay3");
\r
156 h1->GetYaxis()->SetTitle("Counts");
\r
157 fAliITSQADataMakerRec->Add2RawsList((new TH1D(*h1)),1+fGenOffset);
\r
161 if(fLDC==0 || fLDC==3 || fLDC==4){
\r
162 TH1D *h2 = new TH1D("LadPatternL4","Ladder pattern L4",22,0.5,22.5);
\r
163 h2->GetXaxis()->SetTitle("Ladder Number on Lay4");
\r
164 h2->GetYaxis()->SetTitle("Counts");
\r
165 fAliITSQADataMakerRec->Add2RawsList((new TH1D(*h2)),2+fGenOffset);
\r
169 if(fLDC==0 || fLDC==1 || fLDC==2){
\r
170 for(Int_t i=1; i<=fgkLADDonLAY3; i++) {
\r
171 sprintf(hname0,"ModPattern_L3_%d",i);
\r
172 TH1D *h3 = new TH1D(hname0,hname0,6,0.5,6.5);
\r
173 h3->GetXaxis()->SetTitle("Module Number");
\r
174 h3->GetYaxis()->SetTitle("Counts");
\r
175 fAliITSQADataMakerRec->Add2RawsList((new TH1D(*h3)),i-1+3+fGenOffset);
\r
180 if(fLDC==0 || fLDC==3 || fLDC==4){
\r
181 for(Int_t i=1; i<=fgkLADDonLAY4; i++) {
\r
182 sprintf(hname0,"ModPattern_L4_%d",i);
\r
183 TH1D *h4 = new TH1D(hname0,hname0,8,0.5,8.5);
\r
184 h4->GetXaxis()->SetTitle("Module Number");
\r
185 h4->GetYaxis()->SetTitle("Counts");
\r
186 fAliITSQADataMakerRec->Add2RawsList((new TH1D(*h4)),i-1+17+fGenOffset);
\r
192 Int_t indexlast1 = 0;
\r
193 Int_t indexlast2 = 0;
\r
199 indexlast1 = fSDDhTask;
\r
202 for(Int_t i=0; i<3; i++) hname[i]= new char[50];
\r
203 for(Int_t moduleSDD =0; moduleSDD<fgknSDDmodules; moduleSDD++){
\r
204 for(Int_t iside=0;iside<fgknSide;iside++){
\r
205 AliITSgeomTGeo::GetModuleId(moduleSDD+fgkmodoffset, lay, lad, det);
\r
206 sprintf(hname[0],"chargeMapFSE_L%d_%d_%d_%d",lay,lad,det,iside);
\r
207 sprintf(hname[1],"ChargeMapForSingleEvent_L%d_%d_%d_%d",lay,lad,det,iside);
\r
208 sprintf(hname[2],"hmonoDMap_L%d_%d_%d_%d",lay,lad,det,iside);
\r
209 TProfile2D *fModuleChargeMapFSE = new TProfile2D(hname[0],hname[1],256/fTimeBinSize,-0.5,255.5,256,-0.5,255.5);
\r
210 fModuleChargeMapFSE->GetXaxis()->SetTitle("Time Bin");
\r
211 fModuleChargeMapFSE->GetYaxis()->SetTitle("Anode");
\r
212 fAliITSQADataMakerRec->Add2RawsList((new TProfile2D(*fModuleChargeMapFSE)),indexlast1 + index1 + fGenOffset);
\r
213 delete fModuleChargeMapFSE;
\r
217 indexlast2 = indexlast1 + index1;
\r
221 for(Int_t moduleSDD =0; moduleSDD<fgknSDDmodules; moduleSDD++){
\r
222 for(Int_t iside=0;iside<fgknSide;iside++){
\r
223 AliITSgeomTGeo::GetModuleId(moduleSDD+fgkmodoffset, lay, lad, det);
\r
224 sprintf(hname[0],"chargeMap_L%d_%d_%d_%d",lay,lad,det,iside);
\r
225 sprintf(hname[1],"ChargeMap_L%d_%d_%d_%d",lay,lad,det,iside);
\r
226 TProfile2D *fModuleChargeMap = new TProfile2D(hname[0],hname[1],256/fTimeBinSize,-0.5,255.5,256,-0.5,255.5);
\r
227 fModuleChargeMap->GetXaxis()->SetTitle("Time Bin");
\r
228 fModuleChargeMap->GetYaxis()->SetTitle("Anode");
\r
229 fAliITSQADataMakerRec->Add2RawsList((new TProfile2D(*fModuleChargeMap)),indexlast1 + index1 + fGenOffset);
\r
230 delete fModuleChargeMap;
\r
234 indexlast2 = indexlast1 + index1;
\r
240 AliDebug(1,Form("%d SDD Raws histograms booked\n",fSDDhTask));
\r
244 //____________________________________________________________________________
\r
245 void AliITSQASDDDataMakerRec::MakeRaws(AliRawReader* rawReader)
\r
247 // Fill QA for RAW - SDD -
\r
248 if(!fDDLModuleMap){
\r
249 AliError("SDD DDL module map not available - skipping SDD QA");
\r
252 if(rawReader->GetType() != 7) return; // skips non physical triggers
\r
253 AliDebug(1,"entering MakeRaws\n");
\r
254 rawReader->SelectEquipment(17,fgkeqOffset,fgkeqOffset + fgkDDLidRange);
\r
256 rawReader->Reset();
\r
257 AliITSRawStreamSDD s(rawReader);
\r
258 s.SetDDLModuleMap(fDDLModuleMap);
\r
259 Int_t lay, lad, det;
\r
263 for(Int_t moduleSDD =0; moduleSDD<fgknSDDmodules; moduleSDD++){
\r
264 for(Int_t iside=0;iside<fgknSide;iside++) {
\r
265 if(fSDDhTask > 39 + index) fAliITSQADataMakerRec->GetRawsData(39 + index +fGenOffset)->Reset();
\r
273 Int_t isddmod = -1;
\r
274 Int_t coord1, coord2, signal, moduleSDD, ioffset, iorder, activeModule, index1;
\r
276 ildcID = rawReader->GetLDCId();
\r
277 iddl = rawReader->GetDDLID() - fgkDDLIDshift;
\r
278 isddmod = s.GetModuleNumber(iddl,s.GetCarlosId());
\r
280 AliDebug(1,Form("Found module with iddl: %d, s.GetCarlosId: %d \n",iddl,s.GetCarlosId() ));
\r
283 if(s.IsCompletedModule()) {
\r
284 AliDebug(1,Form("IsCompletedModule == KTRUE\n"));
\r
288 coord1 = s.GetCoord1();
\r
289 coord2 = s.GetCoord2();
\r
290 signal = s.GetSignal();
\r
292 moduleSDD = isddmod - fgkmodoffset;
\r
293 if(moduleSDD < 0 || moduleSDD>fgknSDDmodules) {
\r
294 AliDebug(1,Form( "Module SDD = %d, resetting it to 1 \n",moduleSDD));
\r
297 fAliITSQADataMakerRec->GetRawsData(0 +fGenOffset)->Fill(moduleSDD);
\r
299 AliITSgeomTGeo::GetModuleId(isddmod, lay, lad, det);
\r
306 fAliITSQADataMakerRec->GetRawsData(iorder +fGenOffset)->Fill(lad);
\r
307 fAliITSQADataMakerRec->GetRawsData(ioffset+lad-1 +fGenOffset)->Fill(det); //-1 because ladder# starts from 1
\r
309 Short_t iside = s.GetChannel();
\r
310 activeModule = moduleSDD;
\r
311 index1 = activeModule * 2 + iside;
\r
314 AliDebug(1,Form("Wrong index number %d - patched to 0\n",index1));
\r
319 if(fSDDhTask > 39 + index1) {
\r
320 ((TProfile2D *)(fAliITSQADataMakerRec->GetRawsData(39 + index1 +fGenOffset)))->Fill(coord2, coord1, signal);
\r
321 ((TProfile2D *)(fAliITSQADataMakerRec->GetRawsData(39 + index1 + 260*2 +fGenOffset)))->Fill(coord2, coord1, signal);
\r
325 if(!(cnt%10000)) AliDebug(1,Form(" %d raw digits read",cnt));
\r
327 AliDebug(1,Form("Event completed, %d raw digits read",cnt));
\r
330 //____________________________________________________________________________
\r
331 void AliITSQASDDDataMakerRec::InitRecPoints()
\r
333 // Initialization for RECPOINTS - SDD -
\r
334 fGenOffset = (fAliITSQADataMakerRec->fRecPointsQAList)->GetEntries();
\r
336 TH1F *h0 = new TH1F("Lay3TotCh","Layer 3 total charge",1000,-0.5, 499.5);
\r
337 h0->GetXaxis()->SetTitle("ADC value");
\r
338 h0->GetYaxis()->SetTitle("Entries");
\r
339 fAliITSQADataMakerRec->Add2RecPointsList((new TH1F(*h0)), 0 +fGenOffset);
\r
343 TH1F *h1 = new TH1F("Lay4TotCh","Layer 4 total charge",1000,-0.5, 499.5);
\r
344 h1->GetXaxis()->SetTitle("ADC value");
\r
345 h1->GetYaxis()->SetTitle("Entries");
\r
346 fAliITSQADataMakerRec->Add2RecPointsList((new TH1F(*h1)), 1 +fGenOffset);
\r
352 for(Int_t i=1; i<=3; i++){
\r
353 sprintf(hisnam,"Charge_L3_Strip%d",i);
\r
354 TH1F *h2 = new TH1F(hisnam,hisnam,1000,-0.5, 499.5);
\r
355 fAliITSQADataMakerRec->Add2RecPointsList((new TH1F(*h2)),i+1 +fGenOffset);
\r
360 for(Int_t i=1; i<=4; i++){
\r
361 sprintf(hisnam,"Charge_L4_Strip%d",i);
\r
362 TH1F *h3 = new TH1F(hisnam,hisnam,1000,-0.5, 499.5);
\r
363 fAliITSQADataMakerRec->Add2RecPointsList((new TH1F(*h3)),i+4 +fGenOffset);
\r
368 TH1F *h4 = new TH1F("ModPatternRP","Modules pattern RP",fgknSDDmodules,239.5,499.5);
\r
369 h4->GetXaxis()->SetTitle("Module number");
\r
370 h4->GetYaxis()->SetTitle("Entries");
\r
371 fAliITSQADataMakerRec->Add2RecPointsList((new TH1F(*h4)),9 +fGenOffset);
\r
374 TH1F *h5 = new TH1F("ModPatternL3 RP","Ladder pattern L3 RP",14,0.5,14.5);
\r
375 h5->GetXaxis()->SetTitle("Ladder #, Layer 3");
\r
376 h5->GetYaxis()->SetTitle("Entries");
\r
377 fAliITSQADataMakerRec->Add2RecPointsList((new TH1F(*h5)),10 +fGenOffset);
\r
380 TH1F *h6 = new TH1F("ModPatternL4 RP","Ladder pattern L4 RP",22,0.5,22.5);
\r
381 h6->GetXaxis()->SetTitle("Ladder #, Layer 4");
\r
382 h6->GetYaxis()->SetTitle("Entries");
\r
383 fAliITSQADataMakerRec->Add2RecPointsList((new TH1F(*h6)),11 +fGenOffset);
\r
386 TH2F *h7 = new TH2F("Local Coord Distrib","Local Coord Distrib",1000,-4,4,1000,-4,4);
\r
387 h7->GetXaxis()->SetTitle("X local coord, drift, cm");
\r
388 h7->GetYaxis()->SetTitle("Z local coord, anode, cm");
\r
389 fAliITSQADataMakerRec->Add2RecPointsList((new TH2F(*h7)),12 +fGenOffset);
\r
392 TH2F *h8 = new TH2F("Global Coord Distrib","Global Coord Distrib",6000,-30,30,6000,-30,30);
\r
393 h8->GetYaxis()->SetTitle("Y glob coord, cm");
\r
394 h8->GetXaxis()->SetTitle("X glob coord, cm");
\r
395 fAliITSQADataMakerRec->Add2RecPointsList((new TH2F(*h8)),13 +fGenOffset);
\r
399 for(Int_t iLay=0; iLay<=1; iLay++){
\r
400 sprintf(hisnam,"hr_Layer%d",iLay+3);
\r
401 TH1F *h9 = new TH1F(hisnam,hisnam,100,10.,30.);
\r
402 h9->GetXaxis()->SetTitle("r (cm)");
\r
403 h9->GetXaxis()->CenterTitle();
\r
404 h9->GetYaxis()->SetTitle("Entries");
\r
405 fAliITSQADataMakerRec->Add2RecPointsList((new TH1F(*h9)),iLay+14 +fGenOffset);
\r
410 for(Int_t iLay=0; iLay<=1; iLay++){
\r
411 sprintf(hisnam,"hphi_Layer%d",iLay+3);
\r
412 TH1F *h10 = new TH1F(hisnam,hisnam,100,-TMath::Pi(),TMath::Pi());
\r
413 h10->GetXaxis()->SetTitle("#varphi (rad)");
\r
414 h10->GetXaxis()->CenterTitle();
\r
415 h10->GetYaxis()->SetTitle("Entries");
\r
416 fAliITSQADataMakerRec->Add2RecPointsList((new TH1F(*h10)),iLay+16 +fGenOffset);
\r
421 AliDebug(1,Form("%d SDD Recs histograms booked\n",fSDDhTask));
\r
424 //____________________________________________________________________________
\r
425 void AliITSQASDDDataMakerRec::MakeRecPoints(TTree * clustersTree)
\r
427 // Fill QA for RecPoints - SDD -
\r
428 Int_t lay, lad, det;
\r
429 TBranch *branchRecP = clustersTree->GetBranch("ITSRecPoints");
\r
430 if (!branchRecP) {
\r
431 AliError("can't get the branch with the ITS clusters !");
\r
434 static TClonesArray statRecpoints("AliITSRecPoint") ;
\r
435 TClonesArray *recpoints = &statRecpoints;
\r
436 branchRecP->SetAddress(&recpoints);
\r
437 Int_t npoints = 0;
\r
438 Float_t cluglo[3]={0.,0.,0.};
\r
439 for(Int_t module=0; module<clustersTree->GetEntries();module++){
\r
440 branchRecP->GetEvent(module);
\r
441 npoints += recpoints->GetEntries();
\r
442 AliITSgeomTGeo::GetModuleId(module, lay, lad, det);
\r
443 //printf("modnumb %d, lay %d, lad %d, det %d \n",module, lay, lad, det);
\r
445 for(Int_t j=0;j<recpoints->GetEntries();j++){
\r
446 AliITSRecPoint *recp = (AliITSRecPoint*)recpoints->At(j);
\r
447 fAliITSQADataMakerRec->GetRecPointsData(9 +fGenOffset)->Fill(module);
\r
448 recp->GetGlobalXYZ(cluglo);
\r
449 Float_t rad=TMath::Sqrt(cluglo[0]*cluglo[0]+cluglo[1]*cluglo[1]);
\r
450 Float_t phi=TMath::ATan2(cluglo[1],cluglo[0]);
\r
451 if(recp->GetLayer() ==2) {
\r
452 fAliITSQADataMakerRec->GetRecPointsData(0 +fGenOffset)->Fill(recp->GetQ()) ;
\r
453 fAliITSQADataMakerRec->GetRecPointsData(10 +fGenOffset)->Fill(lad);
\r
454 fAliITSQADataMakerRec->GetRecPointsData(14 +fGenOffset)->Fill(rad);
\r
455 fAliITSQADataMakerRec->GetRecPointsData(16 +fGenOffset)->Fill(phi);
\r
456 fAliITSQADataMakerRec->GetRecPointsData(9 +fGenOffset)->Fill(module);
\r
457 fAliITSQADataMakerRec->GetRecPointsData(12 +fGenOffset)->Fill(recp->GetDetLocalX(),recp->GetDetLocalZ());
\r
458 fAliITSQADataMakerRec->GetRecPointsData(13 +fGenOffset)->Fill(cluglo[0],cluglo[1]);
\r
460 else if(recp->GetLayer() ==3) {
\r
461 fAliITSQADataMakerRec->GetRecPointsData(1 +fGenOffset)->Fill(recp->GetQ()) ;
\r
462 fAliITSQADataMakerRec->GetRecPointsData(11 +fGenOffset)->Fill(lad);
\r
463 fAliITSQADataMakerRec->GetRecPointsData(15 +fGenOffset)->Fill(rad);
\r
464 fAliITSQADataMakerRec->GetRecPointsData(17 +fGenOffset)->Fill(phi);
\r
465 fAliITSQADataMakerRec->GetRecPointsData(9 +fGenOffset)->Fill(module);
\r
466 fAliITSQADataMakerRec->GetRecPointsData(12 +fGenOffset)->Fill(recp->GetDetLocalX(),recp->GetDetLocalZ());
\r
467 fAliITSQADataMakerRec->GetRecPointsData(13 +fGenOffset)->Fill(cluglo[0],cluglo[1]);
\r
471 statRecpoints.Clear();
\r