1 #if !defined(__CINT__) || defined(__MAKECINT__)
10 #include "Riostream.h"
12 #include "AliITSgeomTGeo.h"
13 #include "AliITSGainSSDv2.h"
14 #include "AliITSBadChannelsSSDv2.h"
15 #include "AliITSNoiseSSDv2.h"
16 #include "AliITSPedestalSSDv2.h"
17 #include "AliITSGainSSD.h"
18 #include "AliITSBadChannelsSSD.h"
19 #include "AliITSNoiseSSD.h"
20 #include "AliITSPedestalSSD.h"
21 #include "AliCDBManager.h"
22 #include "AliCDBEntry.h"
29 //====================================================================//
30 void Noise(AliCDBManager * man);
31 void Pedestal(AliCDBManager * man);
32 void BadChannelMap(AliCDBManager * man);
33 void GainCalibration(AliCDBManager * man);
34 void ReadOldSSDPedestal(TObjArray *array, AliITSPedestalSSDv2 *pedestalSSD);
35 void ReadOldSSDNoise(TObjArray *array, AliITSNoiseSSDv2 *noiseSSD);
36 void ReadOldSSDBadChannels(TObjArray *array, AliITSBadChannelsSSDv2 *badChannelsSSD);
37 void ReadOldSSDGain(TObjArray *array, AliITSGainSSDv2 *gainSSD);
38 void drawNoiseDistributions(Int_t runNumber);
39 //====================================================================//
41 //_____________________________________________________________________//
42 void readSSDOCDBEntry(const char* type = "alien", Int_t runNumber = 0) {
43 //This macro allows to visualize the bad channels in the OCDB
44 //The type can be either "local" or "alien" (where the OCDB file comes from)
45 //The run nmber is the pedestal one
46 gStyle->SetPalette(1,0);
49 AliCDBManager * man = AliCDBManager::Instance();
51 if(gType == "alien") {
52 man->SetDefaultStorage("alien://folder=/alice/data/2009/Reference/");
54 else if(gType == "local")
55 man->SetDefaultStorage("local://$ALICE_ROOT/OCDB");
57 cout<<"Allowed types: local or alien!!!"<<endl;
61 man->SetRun(runNumber);
66 //GainCalibration(man);
69 //_____________________________________________________________________//
70 void drawNoiseDistributions(Int_t runNumber) {
71 //Draws the noise distributions for each side and layer
72 TString filename = "noiseDistributionsSSD."; filename += runNumber;
74 TFile *f = TFile::Open(filename.Data());
75 TH1F *gHistNoisePSideLayer5 = dynamic_cast<TH1F *>(f->Get("gHistNoisePSideLayer5"));
76 TH1F *gHistNoiseNSideLayer5 = dynamic_cast<TH1F *>(f->Get("gHistNoiseNSideLayer5"));
77 TH1F *gHistNoisePSideLayer6 = dynamic_cast<TH1F *>(f->Get("gHistNoisePSideLayer6"));
78 TH1F *gHistNoiseNSideLayer6 = dynamic_cast<TH1F *>(f->Get("gHistNoiseNSideLayer6"));
80 TCanvas *c1 = new TCanvas("c1","Noise distribution (P-side, Layer 5)",
82 c1->SetFillColor(10); c1->SetHighLightColor(10); c1->SetLogy();
83 gHistNoisePSideLayer5->SetStats(kFALSE); gHistNoisePSideLayer5->Draw();
85 TCanvas *c2 = new TCanvas("c2","Noise distribution (N-side, Layer 5)",
87 c2->SetFillColor(10); c2->SetHighLightColor(10); c2->SetLogy();
88 gHistNoiseNSideLayer5->SetStats(kFALSE); gHistNoiseNSideLayer5->Draw();
90 TCanvas *c3 = new TCanvas("c3","Noise distribution (P-side, Layer 6)",
92 c3->SetFillColor(10); c3->SetHighLightColor(10); c3->SetLogy();
93 gHistNoisePSideLayer6->SetStats(kFALSE); gHistNoisePSideLayer6->Draw();
95 TCanvas *c4 = new TCanvas("c4","Noise distribution (N-side, Layer 6)",
97 c4->SetFillColor(10); c4->SetHighLightColor(10); c4->SetLogy();
98 gHistNoiseNSideLayer6->SetStats(kFALSE); gHistNoiseNSideLayer6->Draw();
101 //_____________________________________________________________________//
102 void drawPedestalDistributions(Int_t runNumber) {
103 //Draws the pedestal distributions for each side and layer
104 TString filename = "pedestalDistributionsSSD."; filename += runNumber;
106 TFile *f = TFile::Open(filename.Data());
107 TH1F *gHistPedestalPSideLayer5 = dynamic_cast<TH1F *>(f->Get("gHistPedestalPSideLayer5"));
108 TH1F *gHistPedestalNSideLayer5 = dynamic_cast<TH1F *>(f->Get("gHistPedestalNSideLayer5"));
109 TH1F *gHistPedestalPSideLayer6 = dynamic_cast<TH1F *>(f->Get("gHistPedestalPSideLayer6"));
110 TH1F *gHistPedestalNSideLayer6 = dynamic_cast<TH1F *>(f->Get("gHistPedestalNSideLayer6"));
112 TCanvas *c1 = new TCanvas("c1","Pedestal distribution (P-side, Layer 5)",
114 c1->SetFillColor(10); c1->SetHighLightColor(10); c1->SetLogy();
115 gHistPedestalPSideLayer5->SetStats(kFALSE); gHistPedestalPSideLayer5->Draw();
117 TCanvas *c2 = new TCanvas("c2","Pedestal distribution (N-side, Layer 5)",
119 c2->SetFillColor(10); c2->SetHighLightColor(10); c2->SetLogy();
120 gHistPedestalNSideLayer5->SetStats(kFALSE); gHistPedestalNSideLayer5->Draw();
122 TCanvas *c3 = new TCanvas("c3","Pedestal distribution (P-side, Layer 6)",
124 c3->SetFillColor(10); c3->SetHighLightColor(10); c3->SetLogy();
125 gHistPedestalPSideLayer6->SetStats(kFALSE); gHistPedestalPSideLayer6->Draw();
127 TCanvas *c4 = new TCanvas("c4","Pedestal distribution (N-side, Layer 6)",
129 c4->SetFillColor(10); c4->SetHighLightColor(10); c4->SetLogy();
130 gHistPedestalNSideLayer6->SetStats(kFALSE); gHistPedestalNSideLayer6->Draw();
133 //_____________________________________________________________________//
134 void Pedestal(AliCDBManager * man) {
135 //Reads the noise OCDB file
136 const Int_t fgkSSDMODULES = 1698;
137 const Int_t fgkSSDSTRIPSPERMODULE = 1536;
138 static const Int_t fgkDefaultNStripsSSD = 768;
140 Int_t runNumber = man->GetRun();
142 //pedestal histograms
143 TH1F *gHistPedestalPSideLayer5 = new TH1F("gHistPedestalPSideLayer5",
144 "Pedestal values (P-side, Layer5); ADC counts; Entries;",
146 TH1F *gHistPedestalNSideLayer5 = new TH1F("gHistPedestalNSideLayer5",
147 "Pedestal values (N-side, Layer5); ADC counts; Entries;",
149 TH1F *gHistPedestalPSideLayer6 = new TH1F("gHistPedestalPSideLayer6",
150 "Pedestal values (P-side, Layer6); ADC counts; Entries;",
152 TH1F *gHistPedestalNSideLayer6 = new TH1F("gHistPedestalNSideLayer6",
153 "Pedestal values (N-side, Layer6); ADC counts; Entries;",
156 Int_t fLayer = 0,fLadder = 0, fModule = 0;
158 AliITSPedestalSSDv2 *pedestalSSD = new AliITSPedestalSSDv2();
159 AliCDBEntry *entryPedestalSSD = man->Get("ITS/Ref/PedestalSSD");
160 TObject *empty = (TObject *)entryPedestalSSD->GetObject();
161 TString objectname = empty->GetName();
162 if(objectname=="TObjArray") {
163 TObjArray *pedestalSSDOld = (TObjArray *)entryPedestalSSD->GetObject();
164 ReadOldSSDPedestal(pedestalSSDOld, pedestalSSD);
166 else if(objectname=="AliITSPedestalSSDv2") {
167 cout<<"Reading the new format of the calibration file..."<<endl;
168 pedestalSSD = (AliITSPedestalSSDv2 *)entryPedestalSSD->GetObject();
171 Double_t pedestal = 0.0;
172 for (Int_t i = 0; i < fgkSSDMODULES; i++) {
173 AliITSgeomTGeo::GetModuleId(i+500,fLayer,fLadder,fModule);
174 //cout<<"Pedestal for module: "<<i+500<<" - Layer: "<<fLayer<<endl;
175 for(Int_t j = 0; j < fgkDefaultNStripsSSD; j++) {
176 pedestal = pedestalSSD->GetPedestalP(i,j);
178 gHistPedestalPSideLayer5->Fill(pedestal);
180 gHistPedestalPSideLayer6->Fill(pedestal);
182 pedestal = pedestalSSD->GetPedestalN(i,j);
184 gHistPedestalNSideLayer5->Fill(pedestal);
186 gHistPedestalNSideLayer6->Fill(pedestal);
190 TString output = "pedestalDistributionsSSD."; output += runNumber;
192 TFile *f = TFile::Open(output.Data(),"recreate");
193 gHistPedestalPSideLayer5->Write();
194 gHistPedestalNSideLayer5->Write();
195 gHistPedestalPSideLayer6->Write();
196 gHistPedestalNSideLayer6->Write();
200 //_____________________________________________________________________//
201 void Noise(AliCDBManager * man) {
202 //Reads the noise OCDB file
203 const Int_t fgkSSDMODULES = 1698;
204 const Int_t fgkSSDSTRIPSPERMODULE = 1536;
205 static const Int_t fgkDefaultNStripsSSD = 768;
207 Int_t runNumber = man->GetRun();
210 TH1F *gHistNoisePSideLayer5 = new TH1F("gHistNoisePSideLayer5",
211 "Noise values (P-side, Layer5); ADC counts; Entries;",
213 TH1F *gHistNoiseNSideLayer5 = new TH1F("gHistNoiseNSideLayer5",
214 "Noise values (N-side, Layer5); ADC counts; Entries;",
216 TH1F *gHistNoisePSideLayer6 = new TH1F("gHistNoisePSideLayer6",
217 "Noise values (P-side, Layer6); ADC counts; Entries;",
219 TH1F *gHistNoiseNSideLayer6 = new TH1F("gHistNoiseNSideLayer6",
220 "Noise values (N-side, Layer6); ADC counts; Entries;",
223 Int_t fLayer = 0,fLadder = 0, fModule = 0;
224 Int_t fHistCounter = 0;
226 TObjArray *array = new TObjArray();
227 TH1D *hNoiseModule[fgkSSDMODULES];
228 for(Int_t i = 500; i < fgkSSDMODULES + 500; i++) {
229 AliITSgeomTGeo::GetModuleId(i,fLayer,fLadder,fModule);
230 fTitle = "SSD_Noise_Layer"; fTitle += fLayer;
231 fTitle += "_Ladder"; fTitle += fLadder;
232 fTitle += "_Module"; fTitle += fModule;
234 hNoiseModule[fHistCounter] = new TH1D(fTitle.Data(),fTitle.Data(),1540,0,1540);
235 hNoiseModule[fHistCounter]->GetXaxis()->SetTitleColor(1);
236 hNoiseModule[fHistCounter]->GetXaxis()->SetTitle("Strip number");
237 hNoiseModule[fHistCounter]->GetYaxis()->SetTitle("Noise");
238 array->AddLast(hNoiseModule[fHistCounter]);
242 AliITSNoiseSSDv2 *noiseSSD = new AliITSNoiseSSDv2();
243 AliCDBEntry *entryNoiseSSD = man->Get("ITS/Calib/NoiseSSD");
244 TObject *empty = (TObject *)entryNoiseSSD->GetObject();
245 TString objectname = empty->GetName();
246 if(objectname=="TObjArray") {
247 TObjArray *noiseSSDOld = (TObjArray *)entryNoiseSSD->GetObject();
248 ReadOldSSDNoise(noiseSSDOld, noiseSSD);
250 else if(objectname=="AliITSNoiseSSDv2") {
251 cout<<"Reading the new format of the calibration file..."<<endl;
252 noiseSSD = (AliITSNoiseSSDv2 *)entryNoiseSSD->GetObject();
255 Double_t noise = 0.0;
256 for (Int_t i = 0; i < fgkSSDMODULES; i++) {
257 AliITSgeomTGeo::GetModuleId(i+500,fLayer,fLadder,fModule);
258 //cout<<"Noise for module: "<<i+500<<" - Layer: "<<fLayer<<endl;
259 for(Int_t j = 0; j < fgkDefaultNStripsSSD; j++) {
260 noise = noiseSSD->GetNoiseP(i,j);
261 hNoiseModule[i]->SetBinContent(j+1,noise);
263 gHistNoisePSideLayer5->Fill(noise);
265 gHistNoisePSideLayer6->Fill(noise);
267 noise = noiseSSD->GetNoiseN(i,j);
268 hNoiseModule[i]->SetBinContent(fgkSSDSTRIPSPERMODULE-j,noise);
270 gHistNoiseNSideLayer5->Fill(noise);
272 gHistNoiseNSideLayer6->Fill(noise);
276 TString output1 = "noiseSSD."; output1 += runNumber; output1 += ".root";
277 TFile *f1 = TFile::Open(output1.Data(),"recreate");
281 TString output2 = "noiseDistributionsSSD."; output2 += runNumber;
283 TFile *f2 = TFile::Open(output2.Data(),"recreate");
284 gHistNoisePSideLayer5->Write();
285 gHistNoiseNSideLayer5->Write();
286 gHistNoisePSideLayer6->Write();
287 gHistNoiseNSideLayer6->Write();
291 //_____________________________________________________________________//
292 void BadChannelMap(AliCDBManager * man) {
293 const Int_t fgkSSDMODULES = 1698;
294 static const Int_t fgkDefaultNStripsSSD = 768;
296 //_____________________________________________________________________________//
297 TH2F *fHistPSideBadChannelMapLayer5 = new TH2F("fHistPSideBadChannelMapLayer5",
298 "Layer 5;N_{module};N_{ladder}",
301 fHistPSideBadChannelMapLayer5->GetXaxis()->SetTitleColor(1);
302 fHistPSideBadChannelMapLayer5->GetZaxis()->SetRangeUser(0.,100.);
303 fHistPSideBadChannelMapLayer5->SetStats(kFALSE);
304 fHistPSideBadChannelMapLayer5->GetYaxis()->SetTitleOffset(1.8);
305 fHistPSideBadChannelMapLayer5->GetXaxis()->SetNdivisions(22);
306 fHistPSideBadChannelMapLayer5->GetYaxis()->SetNdivisions(34);
307 fHistPSideBadChannelMapLayer5->GetXaxis()->SetLabelSize(0.03);
308 fHistPSideBadChannelMapLayer5->GetYaxis()->SetLabelSize(0.03);
309 fHistPSideBadChannelMapLayer5->GetZaxis()->SetTitleOffset(1.6);
310 fHistPSideBadChannelMapLayer5->GetZaxis()->SetTitle("Bad channels (p-side)[%]");
311 TH2F *fHistNSideBadChannelMapLayer5 = new TH2F("fHistNSideBadChannelMapLayer5",
312 "Layer 5;N_{module};N_{ladder}",
315 fHistNSideBadChannelMapLayer5->GetXaxis()->SetTitleColor(1);
316 fHistNSideBadChannelMapLayer5->GetZaxis()->SetRangeUser(0.,100.);
317 fHistNSideBadChannelMapLayer5->SetStats(kFALSE);
318 fHistNSideBadChannelMapLayer5->GetYaxis()->SetTitleOffset(1.8);
319 fHistNSideBadChannelMapLayer5->GetXaxis()->SetNdivisions(22);
320 fHistNSideBadChannelMapLayer5->GetYaxis()->SetNdivisions(34);
321 fHistNSideBadChannelMapLayer5->GetXaxis()->SetLabelSize(0.03);
322 fHistNSideBadChannelMapLayer5->GetYaxis()->SetLabelSize(0.03);
323 fHistNSideBadChannelMapLayer5->GetZaxis()->SetTitleOffset(1.6);
324 fHistNSideBadChannelMapLayer5->GetZaxis()->SetTitle("Bad channels (n-side)[%]");
326 TH2F *fHistPSideBadChannelMapLayer6 = new TH2F("fHistPSideBadChannelMapLayer6",
327 "Layer 6;N_{module};N_{ladder}",
330 fHistPSideBadChannelMapLayer6->GetXaxis()->SetTitleColor(1);
331 fHistPSideBadChannelMapLayer6->GetZaxis()->SetRangeUser(0.,100.);
332 fHistPSideBadChannelMapLayer6->SetStats(kFALSE);
333 fHistPSideBadChannelMapLayer6->GetYaxis()->SetTitleOffset(1.8);
334 fHistPSideBadChannelMapLayer6->GetXaxis()->SetNdivisions(25);
335 fHistPSideBadChannelMapLayer6->GetYaxis()->SetNdivisions(38);
336 fHistPSideBadChannelMapLayer6->GetXaxis()->SetLabelSize(0.03);
337 fHistPSideBadChannelMapLayer6->GetYaxis()->SetLabelSize(0.03);
338 fHistPSideBadChannelMapLayer6->GetZaxis()->SetTitleOffset(1.6);
339 fHistPSideBadChannelMapLayer6->GetZaxis()->SetTitle("Bad channels (p-side)[%]");
340 TH2F *fHistNSideBadChannelMapLayer6 = new TH2F("fHistNSideBadChannelMapLayer6",
341 "Layer 6;N_{module};N_{ladder}",
344 fHistNSideBadChannelMapLayer6->GetXaxis()->SetTitleColor(1);
345 fHistNSideBadChannelMapLayer6->GetZaxis()->SetRangeUser(0.,100.);
346 fHistNSideBadChannelMapLayer6->SetStats(kFALSE);
347 fHistNSideBadChannelMapLayer6->GetYaxis()->SetTitleOffset(1.8);
348 fHistNSideBadChannelMapLayer6->GetXaxis()->SetNdivisions(25);
349 fHistNSideBadChannelMapLayer6->GetYaxis()->SetNdivisions(38);
350 fHistNSideBadChannelMapLayer6->GetXaxis()->SetLabelSize(0.03);
351 fHistNSideBadChannelMapLayer6->GetYaxis()->SetLabelSize(0.03);
352 fHistNSideBadChannelMapLayer6->GetZaxis()->SetTitleOffset(1.6);
353 fHistNSideBadChannelMapLayer6->GetZaxis()->SetTitle("Bad channels (n-side)[%]");
354 //_____________________________________________________________________________//
356 //_____________________________________________________________________________//
357 AliITSBadChannelsSSDv2 *badChannelsSSD = new AliITSBadChannelsSSDv2();
358 AliCDBEntry *entryBadChannelsSSD = man->Get("ITS/Calib/BadChannelsSSD");
359 TObject *empty = (TObject *)entryBadChannelsSSD->GetObject();
360 TString objectname = empty->GetName();
361 if(objectname=="TObjArray") {
362 TObjArray *badChannelsSSDOld = (TObjArray *)entryBadChannelsSSD->GetObject();
363 ReadOldSSDBadChannels(badChannelsSSDOld, badChannelsSSD);
365 else if(objectname=="AliITSBadChannelsSSDv2") {
366 cout<<"Reading the new format of the calibration file..."<<endl;
367 badChannelsSSD = (AliITSBadChannelsSSDv2 *)entryBadChannelsSSD->GetObject();
369 //_____________________________________________________________________________//
371 //_____________________________________________________________________________//
372 Int_t nPSideChannelsTotal = 0, nNSideChannelsTotal = 0;
373 Int_t nBadPSideChannelsTotal = 0, nBadNSideChannelsTotal = 0;
374 Int_t nBadPSideChannels = 0, nBadNSideChannels = 0;
375 Int_t layer = 0, ladder = 0, module = 0;
376 Int_t nPSideChannelsLayer5 = 0, nNSideChannelsLayer5 = 0;
377 Int_t nPSideChannelsLayer6 = 0, nNSideChannelsLayer6 = 0;
378 //_____________________________________________________________________________//
380 for(Int_t i = 0; i < fgkSSDMODULES; i++) {
381 //for(Int_t i = 0; i < 1; i++) {
382 AliITSgeomTGeo::GetModuleId(i+500,layer,ladder,module);
383 nBadPSideChannels = 0, nBadNSideChannels = 0;
384 nPSideChannelsLayer5 = 0, nNSideChannelsLayer5 = 0;
385 nPSideChannelsLayer6 = 0, nNSideChannelsLayer6 = 0;
387 Int_t badChannel = 0;
388 for(Int_t j = 0; j < fgkDefaultNStripsSSD; j++) {
389 badChannel = (Int_t)(badChannelsSSD->GetBadChannelP(i,j));
390 //cout<<"Module: "<<i+500<< " Strip: "<<j<<" - "<<badChannel<<endl;
391 if(badChannel != 0) {
393 nPSideChannelsLayer5 += 1;
395 nPSideChannelsLayer6 += 1;
396 nBadPSideChannels += 1;
398 badChannel = (Int_t)(badChannelsSSD->GetBadChannelN(i,j));
399 //cout<<"Module: "<<i+500<< " Strip: "<<fgkDefaultNStripsSSD+j+1<<" - "<<badChannel<<endl;
400 if(badChannel != 0) {
402 nNSideChannelsLayer5 += 1;
404 nNSideChannelsLayer6 += 1;
405 nBadNSideChannels += 1;
409 if(nPSideChannelsLayer5 > 0)
410 fHistPSideBadChannelMapLayer5->Fill(module,499+ladder,
411 100.*nPSideChannelsLayer5/fgkDefaultNStripsSSD);
412 else fHistPSideBadChannelMapLayer5->Fill(module,499+ladder,0.0001);
413 if(nNSideChannelsLayer5 > 0)
414 fHistNSideBadChannelMapLayer5->Fill(module,499+ladder,
415 100.*nNSideChannelsLayer5/fgkDefaultNStripsSSD);
416 else fHistNSideBadChannelMapLayer5->Fill(module,499+ladder,0.0001);
419 if(nPSideChannelsLayer6 > 0)
420 fHistPSideBadChannelMapLayer6->Fill(module,599+ladder,
421 100.*nPSideChannelsLayer6/fgkDefaultNStripsSSD);
422 else fHistPSideBadChannelMapLayer6->Fill(module,599+ladder,0.0001);
423 if(nNSideChannelsLayer6 > 0)
424 fHistNSideBadChannelMapLayer6->Fill(module,599+ladder,
425 100.*nNSideChannelsLayer6/fgkDefaultNStripsSSD);
426 else fHistNSideBadChannelMapLayer6->Fill(module,599+ladder,0.0001);
429 nBadPSideChannelsTotal += nBadPSideChannels;
430 nBadNSideChannelsTotal += nBadNSideChannels;
431 nPSideChannelsTotal += fgkDefaultNStripsSSD;
432 nNSideChannelsTotal += fgkDefaultNStripsSSD;
435 cout<<"================================="<<endl;
436 cout<<"Bad p-Side channels: "<<100.*nBadPSideChannelsTotal/nPSideChannelsTotal<<endl;
437 cout<<"Bad n-Side channels: "<<100.*nBadNSideChannelsTotal/nNSideChannelsTotal<<endl;
438 cout<<"================================="<<endl;
440 TCanvas *cBadChannel = new TCanvas("cBadChannel",
441 "Bad channel list",0,0,900,900);
442 cBadChannel->SetHighLightColor(10); cBadChannel->SetFillColor(10);
443 cBadChannel->Divide(2,2);
445 cBadChannel->cd(1)->SetBottomMargin(.2);
446 cBadChannel->cd(1)->SetLeftMargin(.15);
447 cBadChannel->cd(1)->SetRightMargin(.2);
448 cBadChannel->cd(1)->SetGridx(); cBadChannel->cd(1)->SetGridy();
449 cBadChannel->cd(1); fHistPSideBadChannelMapLayer5->Draw("colz");
450 cBadChannel->cd(2)->SetBottomMargin(.2);
451 cBadChannel->cd(2)->SetLeftMargin(.15);
452 cBadChannel->cd(2)->SetRightMargin(.2);
453 cBadChannel->cd(2)->SetGridx(); cBadChannel->cd(2)->SetGridy();
454 cBadChannel->cd(2); fHistPSideBadChannelMapLayer6->Draw("colz");
455 cBadChannel->cd(3)->SetBottomMargin(.2);
456 cBadChannel->cd(3)->SetLeftMargin(.15);
457 cBadChannel->cd(3)->SetRightMargin(.2);
458 cBadChannel->cd(3)->SetGridx(); cBadChannel->cd(3)->SetGridy();
459 cBadChannel->cd(3); fHistNSideBadChannelMapLayer5->Draw("colz");
460 cBadChannel->cd(4)->SetBottomMargin(.2);
461 cBadChannel->cd(4)->SetLeftMargin(.15);
462 cBadChannel->cd(4)->SetRightMargin(.2);
463 cBadChannel->cd(4)->SetGridx(); cBadChannel->cd(4)->SetGridy();
464 cBadChannel->cd(4); fHistNSideBadChannelMapLayer6->Draw("colz");
465 cBadChannel->SaveAs("Run-BadChannels.gif");
468 //_____________________________________________________________________//
469 void GainCalibration(AliCDBManager * man) {
470 const Int_t fgkSSDMODULES = 1698;
471 static const Int_t fgkDefaultNStripsSSD = 768;
473 TH2F *fHistGainMapLayer5 = new TH2F("fHistGainMapLayer5",
474 "Layer 5;N_{strip};N_{module}",
477 fHistGainMapLayer5->GetXaxis()->SetTitleColor(1);
478 fHistGainMapLayer5->SetStats(kFALSE);
479 fHistGainMapLayer5->GetYaxis()->SetTitleOffset(1.8);
480 TH2F *fHistGainMapLayer6 = new TH2F("fHistGainMapLayer6",
481 "Layer 6;N_{strip};N_{module}",
484 fHistGainMapLayer6->GetXaxis()->SetTitleColor(1);
485 fHistGainMapLayer6->SetStats(kFALSE);
486 fHistGainMapLayer6->GetYaxis()->SetTitleOffset(1.8);
488 AliITSGainSSDv2 *gainSSD = new AliITSGainSSDv2();
489 AliCDBEntry *entryGainSSD = man->Get("ITS/Calib/GainSSD");
490 TObject *empty = (TObject *)entryGainSSD->GetObject();
491 TString objectname = empty->GetName();
492 if(objectname=="Gain") {
493 TObjArray *gainSSDOld = (TObjArray *)entryGainSSD->GetObject();
494 ReadOldSSDGain(gainSSDOld, gainSSD);
496 else if(objectname=="AliITSGainSSDv2") {
497 cout<<"Reading the new format of the calibration file..."<<endl;
498 gainSSD = (AliITSGainSSDv2 *)entryGainSSD->GetObject();
501 Int_t layer = 0, ladder = 0, module = 0;
503 for (Int_t i = 0; i < fgkSSDMODULES; i++) {
504 //for (Int_t i = 0; i < 1; i++) {
505 AliITSgeomTGeo::GetModuleId(i+500,layer,ladder,module);
506 for(Int_t j = 0; j < fgkDefaultNStripsSSD; j++) {
507 gain = gainSSD->GetGainP(i,j);
508 //cout<<"GainP: "<<gain<<endl;
510 fHistGainMapLayer5->Fill(j,i+500,gain);
512 fHistGainMapLayer6->Fill(j,i+500,gain);
514 gain = gainSSD->GetGainN(i,j);
515 //cout<<"GainN: "<<gain<<endl;
517 fHistGainMapLayer5->Fill(fgkDefaultNStripsSSD+j,i+500,gain);
519 fHistGainMapLayer6->Fill(fgkDefaultNStripsSSD+j,i+500,gain);
523 TCanvas *cGain = new TCanvas("cGain","Gain calibration map",0,300,600,300);
524 cGain->SetHighLightColor(10); cGain->SetFillColor(10); cGain->Divide(2,1);
526 cGain->cd(1)->SetBottomMargin(.2); cGain->cd(1)->SetLeftMargin(.15);
527 cGain->cd(1); fHistGainMapLayer5->Draw("colz");
528 cGain->cd(2)->SetBottomMargin(.2); cGain->cd(2)->SetLeftMargin(.15);
529 cGain->cd(2); fHistGainMapLayer6->Draw("colz");
532 //_____________________________________________________________________//
533 void ReadOldSSDNoise(TObjArray *array,
534 AliITSNoiseSSDv2 *noiseSSD) {
535 const Int_t fgkSSDSTRIPSPERMODULE = 1536;
536 const Int_t fgkSSDPSIDESTRIPSPERMODULE = 768;
538 Int_t fNMod = array->GetEntries();
539 cout<<"Converting old calibration object for noise..."<<endl;
542 Double_t noise = 0.0;
543 for (Int_t iModule = 0; iModule < fNMod; iModule++) {
544 AliITSNoiseSSD *noiseModule = (AliITSNoiseSSD*) (array->At(iModule));
545 for(Int_t iStrip = 0; iStrip < fgkSSDSTRIPSPERMODULE; iStrip++) {
546 noise = (iStrip < fgkSSDPSIDESTRIPSPERMODULE) ? noiseModule->GetNoiseP(iStrip) : noiseModule->GetNoiseN(1535 - iStrip);
547 if(iStrip < fgkSSDPSIDESTRIPSPERMODULE)
548 noiseSSD->AddNoiseP(iModule,iStrip,noise);
549 if(iStrip >= fgkSSDPSIDESTRIPSPERMODULE)
550 noiseSSD->AddNoiseN(iModule,1535 - iStrip,noise);
555 //_____________________________________________________________________//
556 void ReadOldSSDPedestal(TObjArray *array,
557 AliITSPedestalSSDv2 *pedestalSSD) {
558 const Int_t fgkSSDSTRIPSPERMODULE = 1536;
559 const Int_t fgkSSDPSIDESTRIPSPERMODULE = 768;
561 Int_t fNMod = array->GetEntries();
562 cout<<"Converting old calibration object for pedestal..."<<endl;
565 Double_t pedestal = 0.0;
566 for (Int_t iModule = 0; iModule < fNMod; iModule++) {
567 AliITSPedestalSSD *pedestalModule = (AliITSPedestalSSD*) (array->At(iModule));
568 for(Int_t iStrip = 0; iStrip < fgkSSDSTRIPSPERMODULE; iStrip++) {
569 pedestal = (iStrip < fgkSSDPSIDESTRIPSPERMODULE) ? pedestalModule->GetPedestalP(iStrip) : pedestalModule->GetPedestalN(1535 - iStrip);
570 if(iStrip < fgkSSDPSIDESTRIPSPERMODULE)
571 pedestalSSD->AddPedestalP(iModule,iStrip,pedestal);
572 if(iStrip >= fgkSSDPSIDESTRIPSPERMODULE)
573 pedestalSSD->AddPedestalN(iModule,1535 - iStrip,pedestal);
578 //_____________________________________________________________________//
579 void ReadOldSSDBadChannels(TObjArray *array,
580 AliITSBadChannelsSSDv2 *badChannelsSSD) {
581 Int_t fNMod = array->GetEntries();
582 cout<<"Converting old calibration object for bad channels..."<<endl;
583 for (Int_t iModule = 0; iModule < fNMod; iModule++) {
584 //for (Int_t iModule = 0; iModule < 1; iModule++) {
585 AliITSBadChannelsSSD *bad = (AliITSBadChannelsSSD*) (array->At(iModule));
586 TArrayI arrayPSide = bad->GetBadPChannelsList();
587 for(Int_t iPCounter = 0; iPCounter < arrayPSide.GetSize(); iPCounter++)
588 badChannelsSSD->AddBadChannelP(iModule,
590 (Char_t)arrayPSide.At(iPCounter));
592 TArrayI arrayNSide = bad->GetBadNChannelsList();
593 for(Int_t iNCounter = 0; iNCounter < arrayNSide.GetSize(); iNCounter++)
594 badChannelsSSD->AddBadChannelN(iModule,
596 (Char_t)arrayNSide.At(iNCounter));
601 //_____________________________________________________________________//
602 void ReadOldSSDGain(TObjArray *array,
603 AliITSGainSSDv2 *gainSSD) {
604 Int_t fNMod = array->GetEntries();
605 cout<<"Converting old calibration object for gain..."<<endl;
608 for (Int_t iModule = 0; iModule < fNMod; iModule++) {
609 AliITSGainSSD *gainModule = (AliITSGainSSD*) (array->At(iModule));
610 TArrayF arrayPSide = gainModule->GetGainP();
611 for(Int_t iPCounter = 0; iPCounter < arrayPSide.GetSize(); iPCounter++)
612 gainSSD->AddGainP(iModule,
614 arrayPSide.At(iPCounter));
615 TArrayF arrayNSide = gainModule->GetGainN();
616 for(Int_t iNCounter = 0; iNCounter < arrayNSide.GetSize(); iNCounter++)
617 gainSSD->AddGainN(iModule,
619 arrayNSide.At(iNCounter));