initialisation warning fixed
[u/mrichter/AliRoot.git] / PMD / AliPMDCalibGain.cxx
CommitLineData
4fc06b3a 1//////////////////////////////////////////////////////////////////////////////
2//
3// * Copyright(c) 1998-1999, ALICE Experiment at CERN, All rights reserved. *
4// * *
5// * Author: The ALICE Off-line Project. *
6// * Contributors are mentioned in the code where appropriate. *
7// * *
8// * Permission to use, copy, modify and distribute this software and its *
9// * documentation strictly for non-commercial purposes is hereby granted *
10// * without fee, provided that the above copyright notice appears in all *
11// * copies and that both the copyright notice and this permission notice *
12// * appear in the supporting documentation. The authors make no claims *
13// * about the suitability of this software for any purpose. It is *
14// * provided "as is" without express or implied warranty. *
15// **************************************************************************/
16//
17//////////////////////////////////////////////////////////////////////////////
18
19#include "TF1.h"
20#include "TFile.h"
21#include "TObjString.h"
22#include "TROOT.h"
23#include "TClonesArray.h"
24#include "TH1F.h"
25#include "TObjArray.h"
26#include "TTree.h"
7c11ad1f 27#include "TMath.h"
4fc06b3a 28
29// --- Standard library ---
30
31// --- AliRoot header files ---
32#include "AliDAQ.h"
33#include "AliLog.h"
34#include "AliRawReader.h"
35#include "AliPMDRawStream.h"
36#include "AliPMDddldata.h"
37#include "AliPMDCalibGain.h"
38
39ClassImp(AliPMDCalibGain)
40
1da50520 41AliPMDCalibGain::AliPMDCalibGain():
42 TObject(),
43 fpw(NULL)
4fc06b3a 44{
45 // Standard Constructor
7c11ad1f 46 for(Int_t idet = 0; idet < kDet; idet++)
4fc06b3a 47 {
7c11ad1f 48 for(Int_t ismn = 0; ismn < kMaxSMN; ismn++)
4fc06b3a 49 {
7c11ad1f 50 fSMIso[idet][ismn] = 0.;
51 fSMCount[idet][ismn] = 0.;
52 fCountSm[idet][ismn]=0.;
53 fTempnhit[idet][ismn]=0.;
54 fTempnhitSq[idet][ismn]=0.;
55 for(Int_t jrow = 0; jrow < kMaxRow; jrow++)
4fc06b3a 56 {
7c11ad1f 57 for(Int_t kcol = 0; kcol < kMaxCol; kcol++)
4fc06b3a 58 {
7c11ad1f 59 fCellIso[idet][ismn][jrow][kcol] = 0.;
60 fCellCount[idet][ismn][jrow][kcol] = 0.;
61 fNhitCell[idet][ismn][jrow][kcol] = 0.;
62 fPedMeanRMS[idet][ismn][jrow][kcol] = 0.;
63 fHotFlag[idet][ismn][jrow][kcol] = 0.;
64
4fc06b3a 65 }
66 }
67 }
68 }
7c11ad1f 69
4fc06b3a 70}
71// ------------------------------------------------------------------------ //
72AliPMDCalibGain::AliPMDCalibGain(const AliPMDCalibGain &pmdcalibgain):
1da50520 73 TObject(pmdcalibgain),
74 fpw(NULL)
4fc06b3a 75{
7c11ad1f 76 for(Int_t idet = 0; idet < kDet; idet++)
4fc06b3a 77 {
7c11ad1f 78 for(Int_t ismn = 0; ismn < kMaxSMN; ismn++)
4fc06b3a 79 {
7c11ad1f 80 fSMIso[idet][ismn] = pmdcalibgain.fSMIso[idet][ismn] ;
81 fSMCount[idet][ismn] = pmdcalibgain.fSMCount[idet][ismn] ;
82 fCountSm[idet][ismn] = pmdcalibgain.fCountSm[idet][ismn];
83 fTempnhit[idet][ismn] = pmdcalibgain.fTempnhit[idet][ismn];
84 fTempnhitSq[idet][ismn] = pmdcalibgain.fTempnhitSq[idet][ismn];
85 for(Int_t jrow = 0; jrow < kMaxRow; jrow++)
4fc06b3a 86 {
7c11ad1f 87 for(Int_t kcol = 0; kcol < kMaxCol; kcol++)
4fc06b3a 88 {
7c11ad1f 89 fCellIso[idet][ismn][jrow][kcol] = pmdcalibgain.fCellIso[idet][ismn][jrow][kcol];
90 fCellCount[idet][ismn][jrow][kcol] = pmdcalibgain.fCellCount[idet][ismn][jrow][kcol];
91 fNhitCell[idet][ismn][jrow][kcol] = pmdcalibgain.fNhitCell[idet][ismn][jrow][kcol];
92 fPedMeanRMS[idet][ismn][jrow][kcol] = pmdcalibgain.fPedMeanRMS[idet][ismn][jrow][kcol];
93 fHotFlag[idet][ismn][jrow][kcol] = pmdcalibgain.fHotFlag[idet][ismn][jrow][kcol];
94
4fc06b3a 95 }
96 }
97 }
98 }
7c11ad1f 99
4fc06b3a 100}
101// ------------------------------------------------------------------------ //
102AliPMDCalibGain &AliPMDCalibGain::operator=(const AliPMDCalibGain &pmdcalibgain)
103{
7c11ad1f 104 if(this != &pmdcalibgain)
4fc06b3a 105 {
1da50520 106 this->fpw = pmdcalibgain.fpw;
7c11ad1f 107 for(Int_t idet = 0; idet < kDet; idet++)
4fc06b3a 108 {
7c11ad1f 109 for(Int_t ismn = 0; ismn < kMaxSMN; ismn++)
4fc06b3a 110 {
7c11ad1f 111 fSMIso[idet][ismn] = pmdcalibgain.fSMIso[idet][ismn];
112 fSMCount[idet][ismn] = pmdcalibgain.fSMCount[idet][ismn];
113 fCountSm[idet][ismn] = pmdcalibgain.fCountSm[idet][ismn];
114 fTempnhit[idet][ismn] = pmdcalibgain.fTempnhit[idet][ismn];
115 fTempnhitSq[idet][ismn] = pmdcalibgain.fTempnhitSq[idet][ismn];
116 for(Int_t jrow = 0; jrow < kMaxRow;jrow++)
4fc06b3a 117 {
7c11ad1f 118 for(Int_t kcol = 0; kcol < kMaxCol; kcol++)
119 {
120 fCellIso[idet][ismn][jrow][kcol] = pmdcalibgain.fCellIso[idet][ismn][jrow][kcol];
121 fCellCount[idet][ismn][jrow][kcol] = pmdcalibgain.fCellCount[idet][ismn][jrow][kcol];
122 fNhitCell[idet][ismn][jrow][kcol] = pmdcalibgain.fNhitCell[idet][ismn][jrow][kcol]; //za
123 fPedMeanRMS[idet][ismn][jrow][kcol] = pmdcalibgain.fPedMeanRMS[idet][ismn][jrow][kcol];
124 fHotFlag[idet][ismn][jrow][kcol] = pmdcalibgain.fHotFlag[idet][ismn][jrow][kcol];
125
4fc06b3a 126 }
127 }
128 }
129 }
130 }
7c11ad1f 131 return *this;
4fc06b3a 132}
133// ------------------------------------------------------------------------ //
134AliPMDCalibGain::~AliPMDCalibGain()
135{
7c11ad1f 136 // dtor
4fc06b3a 137}
1da50520 138
139// ------------------------------------------------------------------------ //
140
d31b0352 141Int_t AliPMDCalibGain::ExtractPedestal(const Char_t *rootFile)
1da50520 142{
143 // Pedestal extraction from the PMD_PED.root file
144 // To be called once at the beginning
7c11ad1f 145
1da50520 146 Int_t det, sm, row, col;
147 Float_t mean, rms;
7c11ad1f 148
d31b0352 149 TFile *pedfile = new TFile(rootFile);
7c11ad1f 150
1da50520 151 if(!pedfile)
152 {
153 printf("ERROR --- NO PEDESTAL (PMD_PED1.root) FILE IS FOUND --- STOP GAIN DA\n");
154 return -3;
155 }
156
157
158 TTree *ped =(TTree*)pedfile->Get("ped");
159
160 ped->SetBranchAddress("det",&det);
161 ped->SetBranchAddress("sm",&sm);
162 ped->SetBranchAddress("row",&row);
163 ped->SetBranchAddress("col",&col);
164 ped->SetBranchAddress("mean",&mean);
165 ped->SetBranchAddress("rms",&rms);
166
167 Int_t nentries = (Int_t)ped->GetEntries();
168
169 for (Int_t ient = 0; ient < nentries; ient++)
170 {
171 ped->GetEntry(ient);
172 fPedMeanRMS[det][sm][row][col] = mean + 3.*rms;
173 //printf("Mean= %f, RMS= %f, PedMeanRMS=%f\n",mean,rms,fPedMeanRMS[det][sm][row][col]);
174
175 }
176
177 pedfile->Close();
178 delete pedfile;
179 pedfile = 0x0;
180
181 return 1;
182}
7c11ad1f 183//------------------------------------------------------------------------------------------------
184
185Int_t AliPMDCalibGain::ExtractHotChannel(const Char_t *rootFile)
186{
187 // HotChannel extraction from the PMD_HOT.root file
188 // To be called once at the beginning
189
190 Int_t det, sm, row, col;
191 Float_t flag;
192
193 TFile *hotmapfile = new TFile(rootFile);
194
195 if(!hotmapfile)
196 {
9d1486e3 197 printf(" NO HOTCHANNEL MAP FOUND (PMD_HOT.root) FILE IS FOUND \n");
7c11ad1f 198 fHotFlag[kDet][kMaxSMN][kMaxRow][kMaxCol] = 0.;
199 }
200
201
202 TTree *hot =(TTree*)hotmapfile->Get("hot");
203
204 hot->SetBranchAddress("det",&det);
205 hot->SetBranchAddress("sm",&sm);
206 hot->SetBranchAddress("row",&row);
207 hot->SetBranchAddress("col",&col);
208 hot->SetBranchAddress("flag",&flag);
209
210 Int_t nentries = (Int_t)hot->GetEntries();
211
212 for (Int_t ient = 0; ient < nentries; ient++)
213 {
214 hot->GetEntry(ient);
215 fHotFlag[det][sm][row][col] = flag;
216
217 //printf(" HotFlag=%f\n",fHotFlag[det][sm][row][col]);
218 }
219
220 hotmapfile->Close();
221 delete hotmapfile;
222 hotmapfile = 0x0;
223
224 return 1;
225}
226
227
1da50520 228// ------------------------------------------------------------------------ //
229
d31b0352 230void AliPMDCalibGain::ReadTempFile(const Char_t *tempFile)
1da50520 231{
232 // Read the variables from the file
233
d31b0352 234 fpw = fopen(tempFile,"r");
7c11ad1f 235
1da50520 236 Float_t smcount, smiso;
237 Float_t cellcount, celliso;
238
1da50520 239
240 for (Int_t idet = 0; idet < kDet; idet++)
241 {
242 for (Int_t ism = 0; ism < kMaxSMN; ism++)
243 {
244 fscanf(fpw,"%d %d %f %f",&idet,&ism,&smcount,&smiso);
7c11ad1f 245
1da50520 246 fSMCount[idet][ism] = smcount;
247 fSMIso[idet][ism] = smiso;
248 }
249 }
250
251 for (Int_t idet = 0; idet < kDet; idet++)
252 {
253 for (Int_t ism = 0; ism < kMaxSMN; ism++)
254 {
255 for (Int_t irow = 0; irow < kMaxRow; irow++)
256 {
257 for (Int_t icol = 0; icol < kMaxCol; icol++)
258 {
259 fscanf(fpw,"%d %d %d %d %f %f",&idet,&ism,&irow,&icol,
7c11ad1f 260 &cellcount,&celliso);
261
1da50520 262 fCellCount[idet][ism][irow][icol] = cellcount;
263 fCellIso[idet][ism][irow][icol] = celliso;
264 }
265 }
266 }
267 }
7c11ad1f 268
1da50520 269 fclose(fpw);
270
271}
4fc06b3a 272// ------------------------------------------------------------------------ //
d31b0352 273void AliPMDCalibGain::WriteTempFile(const Char_t *tempFile)
1da50520 274{
d31b0352 275 // Write the Temporary file if the required statics is not achieved
1da50520 276
1da50520 277
278 /*
d31b0352 279 Following variables to be written to a file
1da50520 280 fDetIso[idet] ;
281 fSMIso[idet][ismn];
282 fCellIso[idet][ismn][irow][icol];
283
284 fDetCount[idet];
285 fSMCount[idet][ismn];
286 fCellCount[idet][ismn][irow][icol];
287 */
288
d31b0352 289
290 fpw = fopen(tempFile,"w+");
7c11ad1f 291
1da50520 292 for (Int_t idet = 0; idet < kDet; idet++)
293 {
7c11ad1f 294 // fprintf(fpw,"%d %f %f\n",idet,fDetCount[idet],fDetIso[idet]);
1da50520 295 }
7c11ad1f 296
1da50520 297 for (Int_t idet = 0; idet < kDet; idet++)
298 {
299 for (Int_t ism = 0; ism < kMaxSMN; ism++)
300 {
301 fprintf(fpw,"%d %d %f %f\n",idet,ism, fSMCount[idet][ism],fSMIso[idet][ism]);
302 }
303 }
7c11ad1f 304
1da50520 305 for (Int_t idet = 0; idet < kDet; idet++)
306 {
307 for (Int_t ism = 0; ism < kMaxSMN; ism++)
308 {
309 for (Int_t irow = 0; irow < kMaxRow; irow++)
310 {
311 for (Int_t icol = 0; icol < kMaxCol; icol++)
312 {
313 fprintf(fpw,"%d %d %d %d %f %f\n",idet,ism,irow,icol,
314 fCellCount[idet][ism][irow][icol],
315 fCellIso[idet][ism][irow][icol]);
316 }
317 }
318 }
319 }
7c11ad1f 320
1da50520 321 fclose(fpw);
7c11ad1f 322
1da50520 323}
324
325// ------------------------------------------------------------------------ //
326
338c285a 327Bool_t AliPMDCalibGain::ProcessEvent(AliRawReader *rawReader, TObjArray *pmdddlcont)
4fc06b3a 328{
329 // Calculates the ADC of isolated cell
7c11ad1f 330
4fc06b3a 331 const Int_t kDDL = AliDAQ::NumberOfDdls("PMD");
332 const Int_t kCellNeighbour = 6;
4fc06b3a 333 Int_t neibx[6] = {1,0,-1,-1,0,1};
334 Int_t neiby[6] = {0,1,1,0,-1,-1};
335
336 Int_t id1,jd1; //neighbour row/col
337 Int_t isocount; //number of neighbours with 0 signal
338
0ab3a530 339 Float_t d1[kDet][kMaxSMN][kMaxRow][kMaxCol];
7c11ad1f 340
4fc06b3a 341 for(Int_t idet = 0; idet < kDet; idet++)
342 {
343 for(Int_t ismn = 0; ismn < kMaxSMN; ismn++)
344 {
7c11ad1f 345 for(Int_t irow = 0; irow < kMaxRow; irow++)
4fc06b3a 346 {
347 for(Int_t icol = 0; icol < kMaxCol; icol++)
348 {
7c11ad1f 349 d1[idet][ismn][irow][icol] = 0.;
4fc06b3a 350 }
351 }
352 }
353 }
7c11ad1f 354
4fc06b3a 355 AliPMDRawStream rawStream(rawReader);
356
1da50520 357 Int_t iddl = -1;
7c11ad1f 358
722ccc67 359 Int_t numberofDDLs = 0;
360
7c11ad1f 361 while ((iddl = rawStream.DdlData(pmdddlcont)) >=0) {
362 numberofDDLs++;
363
364 Int_t ientries = pmdddlcont->GetEntries();
365
366 for (Int_t ient = 0; ient < ientries; ient++)
4fc06b3a 367 {
7c11ad1f 368 AliPMDddldata *pmdddl = (AliPMDddldata*)pmdddlcont->UncheckedAt(ient);
369
370 Int_t idet = pmdddl->GetDetector();
371 Int_t ismn = pmdddl->GetSMN();
372 Int_t mcm = pmdddl->GetMCM();
373 //Int_t ichno = pmdddl->GetChannel();
374 Int_t irow = pmdddl->GetRow();
375 Int_t icol = pmdddl->GetColumn();
376 Int_t isig = pmdddl->GetSignal();
377
378 // This is the protection not to crash the code
379
380 if(mcm == 0) continue;
381 if (irow < 0 || icol < 0 || irow > 47 || icol > 95) continue;
382
383 // Pedestal subtraction
384 if(fHotFlag[idet][ismn][irow][icol] == 1.0) isig = 0;
385
386 if (isig>0)
387 {
388 d1[idet][ismn][irow][icol] =
389 (Float_t) isig - fPedMeanRMS[idet][ismn][irow][icol];
390
391 //printf("Signal_ped_subtracted=%f, pedestal=%f\n",d1[idet][ismn][irow][icol]),fPedMeanRMS[idet][ismn][irow][icol];
392
393 fNhitCell[idet][ismn][irow][icol]++; //cell hit frequency
394 }
395
1da50520 396
4fc06b3a 397 }
7c11ad1f 398 pmdddlcont->Delete();
399 }
4fc06b3a 400
401 for(Int_t idet=0; idet < kDet; idet++)
7c11ad1f 402 {
4fc06b3a 403 for(Int_t ismn = 0; ismn < kMaxSMN; ismn++)
7c11ad1f 404 {
4fc06b3a 405 for(Int_t irow = 0; irow < kMaxRow; irow++)
7c11ad1f 406 {
4fc06b3a 407 for(Int_t icol = 0; icol < kMaxCol; icol++)
7c11ad1f 408 {
4fc06b3a 409 if(d1[idet][ismn][irow][icol] > 0)
7c11ad1f 410 {
4fc06b3a 411 isocount = 0;
412 for(Int_t ii = 0; ii < kCellNeighbour; ii++)
7c11ad1f 413 {
4fc06b3a 414 id1 = irow + neibx[ii];
415 jd1 = icol + neiby[ii];
d55d21ef 416 if (id1 < 0) id1 = 0;
417 if (id1 > kMaxRow-1) id1 = kMaxRow - 1;
418 if (jd1 < 0) jd1 = 0;
419 if (jd1 > kMaxCol-1) jd1 = kMaxCol - 1;
4fc06b3a 420 if(d1[idet][ismn][id1][jd1] == 0)
7c11ad1f 421 {
4fc06b3a 422 isocount++;
423 if(isocount == kCellNeighbour)
7c11ad1f 424 {
425 //fDetIso[idet] += d1[idet][ismn][irow][icol];
0ab3a530 426 fSMIso[idet][ismn] += d1[idet][ismn][irow][icol];
427 fCellIso[idet][ismn][irow][icol] += d1[idet][ismn][irow][icol];
7c11ad1f 428 //fDetCount[idet]++;
0ab3a530 429 fSMCount[idet][ismn]++;
430 fCellCount[idet][ismn][irow][icol]++;
431
7c11ad1f 432 }
433 }
434 } // neigh cell cond.
435 } // d>0 cond.
436 }
437 }
438 }
439 }
440
441 for(Int_t idet=0; idet < kDet; idet++)
442 {
443 for(Int_t ismn = 0; ismn < kMaxSMN; ismn++)
444 {
445 for(Int_t irow = 0; irow < kMaxRow; irow++)
446 {
447 for(Int_t icol = 0; icol < kMaxCol; icol++)
448 {
449 if(fNhitCell[idet][ismn][irow][icol]>0)
450 {
451 fCountSm[idet][ismn] += 1;
452 fTempnhit[idet][ismn] += fNhitCell[idet][ismn][irow][icol];
453 fTempnhitSq[idet][ismn] += fNhitCell[idet][ismn][irow][icol]
454 *fNhitCell[idet][ismn][irow][icol];
455 }
456 }
457 }
458 }
459 }
460
461 if (numberofDDLs < kDDL) return kFALSE;
722ccc67 462 return kTRUE;
7c11ad1f 463
4fc06b3a 464}
465// ------------------------------------------------------------------------ //
9d1486e3 466void AliPMDCalibGain::Analyse(TTree *gaintree, TTree *meantree)
4fc06b3a 467{
7c11ad1f 468 // Calculates the mean
469 Int_t det, sm, row, col;
470 Float_t gain;
7c11ad1f 471 Float_t cellmean = 0.;
472
9d1486e3 473 Float_t modmean[2][24];
474
475 for (Int_t idet=0; idet < 2; idet++)
476 {
477 for (Int_t ism = 0; ism < 24; ism++)
478 {
479 modmean[idet][ism] = 0.;
480 }
481 }
482
7c11ad1f 483 gaintree->Branch("det",&det,"det/I");
484 gaintree->Branch("sm",&sm,"sm/I");
485 gaintree->Branch("row",&row,"row/I");
486 gaintree->Branch("col",&col,"col/I");
487 gaintree->Branch("gain",&gain,"gain/F");
488
489 for(Int_t idet = 0; idet < kDet; idet++)
490 {
491 for(Int_t ism = 0; ism < kMaxSMN; ism++)
492 {
493 if (fSMCount[idet][ism] > 0)
9d1486e3 494 modmean[idet][ism] = fSMIso[idet][ism]/fSMCount[idet][ism];
7c11ad1f 495 for(Int_t irow = 0; irow < kMaxRow; irow++)
496 {
497 for(Int_t icol = 0; icol < kMaxCol; icol++)
498 {
499 if (fCellCount[idet][ism][irow][icol] > 0.)
500 {
501 cellmean = fCellIso[idet][ism][irow][icol]/fCellCount[idet][ism][irow][icol];
502 }
503 det = idet;
504 sm = ism;
505 row = irow;
506 col = icol;
507 if (cellmean > 0.0 && fCellCount[idet][ism][irow][icol]>0.)
508 {
9d1486e3 509 gain = cellmean/modmean[idet][ism];
7c11ad1f 510 }
511 else
512 {
513 gain = 0.;
514 }
515 //if(fCellCount[idet][ism][irow][icol]>0.) printf("CellCount =%f, gain= %f\n",fCellCount[idet][ism][irow][icol],gain);
516 gaintree->Fill();
517 }
518 }
519 }
520 }
9d1486e3 521
522 Float_t smmean;
523
524 // Writing each module mean value
525 meantree->Branch("det",&det,"det/I");
526 meantree->Branch("sm",&sm,"sm/I");
527 meantree->Branch("smmean",&smmean,"row/F");
7c11ad1f 528
9d1486e3 529 for(Int_t idet = 0; idet < kDet; idet++)
530 {
531 for (Int_t ism = 0; ism < kMaxSMN; ism++)
532 {
533 det = idet;
534 sm = ism;
535 smmean = modmean[idet][ism];
536 meantree->Fill();
537 }
538 }
539
4fc06b3a 540}
541// ------------------------------------------------------------------------ //
9d1486e3 542void AliPMDCalibGain::FindHotCell(TTree *hottree, Float_t xvar)
7c11ad1f 543{
544 // Calculates the mean
545 Int_t det, sm, row, col;
546 Float_t flag;
547 Float_t meannhit;
548 Float_t meanSqnhit;
549 Float_t sigmanhit,nhitcut;
9d1486e3 550
551 //Float_t xvar = 5.;
552
7c11ad1f 553 hottree->Branch("det",&det,"det/I");
554 hottree->Branch("sm",&sm,"sm/I");
555 hottree->Branch("row",&row,"row/I");
556 hottree->Branch("col",&col,"col/I");
557 hottree->Branch("flag",&flag,"flag/F");
558
559 for(Int_t idet = 0; idet < kDet; idet++)
560 {
561 for(Int_t ism = 0; ism < kMaxSMN; ism++)
562 {
563 if (fCountSm[idet][ism]> 0)
564 {
565 meannhit = fTempnhit[idet][ism]/fCountSm[idet][ism];
566 meanSqnhit = fTempnhitSq[idet][ism]/fCountSm[idet][ism];
567 sigmanhit = sqrt(meanSqnhit-(meannhit*meannhit));
9d1486e3 568 nhitcut = meannhit + xvar*sigmanhit;
7c11ad1f 569
570 for(Int_t irow = 0; irow < kMaxRow; irow++)
571 {
572 for(Int_t icol = 0; icol < kMaxCol; icol++)
573 {
574 det = idet;
575 sm = ism;
576 row = irow;
577 col = icol;
578
579 if(fNhitCell[idet][ism][irow][icol] > nhitcut)
580 {
581 flag = 1.0;
582 }
583 else
584 {
585 flag = 0.;
586 }
587 hottree->Fill();
588 }
589
590 }
591 }
592 }
593 }
594}
595