1 /**************************************************************************
2 * Copyright(c) 1998-1999, 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 **************************************************************************/
18 ////////////////////////////////////////////////////////////////////////////
20 // AliTRDCalibraVector //
22 // This class is for the vector method of the TRD calibration. //
25 // R. Bailhache (R.Bailhache@gsi.de) //
27 ////////////////////////////////////////////////////////////////////////////
29 #include <TGraphErrors.h>
31 #include <TObjArray.h>
33 #include <TDirectory.h>
39 #include "AliTRDCalibraVector.h"
40 #include "AliTRDCommonParam.h"
44 ClassImp(AliTRDCalibraVector)
46 //______________________________________________________________________________________
47 AliTRDCalibraVector::AliTRDCalibraVector()
61 // Default constructor
64 for (Int_t idet = 0; idet < 540; idet++){
66 fPHEntries[idet]=new TArrayI();
67 fPHMean[idet]=new TArrayF();
68 fPHSquares[idet]=new TArrayF();
70 fPRFEntries[idet]=new TArrayI();
71 fPRFMean[idet]=new TArrayF();
72 fPRFSquares[idet]=new TArrayF();
75 fCHEntries[idet]=new TArrayI();
79 for(Int_t k = 0; k < 3; k++){
85 //______________________________________________________________________________________
86 AliTRDCalibraVector::AliTRDCalibraVector(const AliTRDCalibraVector &c)
94 ,fNumberBinCharge(c.fNumberBinCharge)
95 ,fNumberBinPRF(c.fNumberBinPRF)
97 ,fPRFRange(c.fPRFRange)
102 for (Int_t idet = 0; idet < 540; idet++){
104 const TArrayI *phEntries = (TArrayI*)c.fPHEntries[idet];
105 const TArrayF *phMean = (TArrayF*)c.fPHMean[idet];
106 const TArrayF *phSquares = (TArrayF*)c.fPHSquares[idet];
108 const TArrayI *prfEntries = (TArrayI*)c.fPRFEntries[idet];
109 const TArrayF *prfMean = (TArrayF*)c.fPRFMean[idet];
110 const TArrayF *prfSquares = (TArrayF*)c.fPRFSquares[idet];
112 const TArrayI *chEntries = (TArrayI*)c.fCHEntries[idet];
114 if ( phEntries != 0x0 ) fPHEntries[idet]=new TArrayI(*phEntries);
115 if ( phMean != 0x0 ) fPHMean[idet]=new TArrayF(*phMean);
116 if ( phSquares != 0x0 ) fPHSquares[idet]=new TArrayF(*phSquares);
118 if ( prfEntries != 0x0 ) fPRFEntries[idet]=new TArrayI(*prfEntries);
119 if ( prfMean != 0x0 ) fPRFMean[idet]=new TArrayF(*prfMean);
120 if ( prfSquares != 0x0 ) fPRFSquares[idet]=new TArrayF(*prfSquares);
123 if ( chEntries != 0x0 ) fCHEntries[idet]=new TArrayI(*chEntries);
127 for(Int_t k = 0; k < 3; k++){
128 fDetCha0[k] = c.fDetCha0[k];
129 fDetCha2[k] = c.fDetCha2[k];
134 //_____________________________________________________________________
135 AliTRDCalibraVector& AliTRDCalibraVector::operator = (const AliTRDCalibraVector &source)
138 // assignment operator
140 if (&source == this) return *this;
141 new (this) AliTRDCalibraVector(source);
145 //____________________________________________________________________________________
146 AliTRDCalibraVector::~AliTRDCalibraVector()
149 // AliTRDCalibraVector destructor
152 if(fPHEntries) delete fPHEntries;
153 if(fPHMean) delete fPHMean;
154 if(fPHSquares) delete fPHSquares;
155 if(fPRFEntries) delete fPRFEntries;
156 if(fPRFMean) delete fPRFMean;
157 if(fPRFSquares) delete fPRFSquares;
158 if(fCHEntries) delete fCHEntries;
161 //_____________________________________________________________________________
162 Int_t AliTRDCalibraVector::SearchBin(Float_t value, Int_t i) const
170 Float_t fbinmax = value;
171 Int_t fNumberOfBin = -1;
178 fNumberOfBin = fNumberBinCharge;
182 fbinmax = TMath::Abs(fPRFRange);
183 fbinmin = -TMath::Abs(fPRFRange);
184 fNumberOfBin = fNumberBinPRF;
192 if ((value >= fbinmax) ||
197 reponse = (Int_t) ((fNumberOfBin*(value-fbinmin)) / (fbinmax-fbinmin));
203 //_____________________________________________________________________________
204 Bool_t AliTRDCalibraVector::UpdateVectorCH(Int_t det, Int_t group, Float_t value)
207 // Fill the vector CH
211 Int_t bin = SearchBin(value,0);
219 if(fDetectorCH != det){
220 fCHEntries[det] = ((TArrayI *)GetCHEntries(det,kTRUE));
223 Int_t entries = ((TArrayI *)fCHEntries[det])->At(group*fNumberBinCharge+bin);
225 Int_t entriesn = entries+1;
226 ((TArrayI *)fCHEntries[det])->AddAt(entriesn,group*fNumberBinCharge+bin);
234 //_____________________________________________________________________________
235 Bool_t AliTRDCalibraVector::UpdateVectorPRF(Int_t det, Int_t group, Float_t x, Float_t y)
238 // Fill the vector PRF
242 Int_t bin = SearchBin(x,2);
249 if(fDetectorPRF != det){
250 fPRFEntries[det] = ((TArrayI *)GetPRFEntries(det,kTRUE));
251 fPRFMean[det] = ((TArrayF *)GetPRFMean(det,kTRUE));
252 fPRFSquares[det] = ((TArrayF *)GetPRFSquares(det,kTRUE));
255 Int_t entries = ((TArrayI *)fPRFEntries[det])->At((Int_t)group*fNumberBinPRF+bin);
256 Float_t mean = ((TArrayI *)fPRFMean[det])->At((Int_t)group*fNumberBinPRF+bin);
257 Float_t square = ((TArrayI *)fPRFSquares[det])->At((Int_t)group*fNumberBinPRF+bin);
259 Int_t entriesn = entries+1;
260 ((TArrayI *)fPRFEntries[det])->AddAt(entriesn,(Int_t)group*fNumberBinPRF+bin);
261 Float_t meann = (mean*((Float_t)entries)+y)/((Float_t)entriesn);
262 ((TArrayF *)fPRFMean[det])->AddAt(meann,(Int_t)group*fNumberBinPRF+bin);
263 Float_t squaren = ((square*((Float_t)entries))+(y*y))/((Float_t)entriesn);
264 ((TArrayF *)fPRFSquares[det])->AddAt(squaren,(Int_t)group*fNumberBinPRF+bin);
272 //_____________________________________________________________________________
273 Bool_t AliTRDCalibraVector::UpdateVectorPH(Int_t det, Int_t group, Int_t time, Float_t value)
276 // Fill the vector PH
288 if(fDetectorPH != det){
289 fPHEntries[det] = ((TArrayI *)GetPHEntries(det,kTRUE));
290 fPHMean[det] = ((TArrayF *)GetPHMean(det,kTRUE));
291 fPHSquares[det] = ((TArrayF *)GetPHSquares(det,kTRUE));
294 Int_t entries = ((TArrayI *)fPHEntries[det])->At(group*fTimeMax+bin);
295 Float_t mean = ((TArrayI *)fPHMean[det])->At(group*fTimeMax+bin);
296 Float_t square = ((TArrayI *)fPHSquares[det])->At(group*fTimeMax+bin);
298 Int_t entriesn = entries+1;
299 ((TArrayI *)fPHEntries[det])->AddAt(entriesn,group*fTimeMax+bin);
300 Float_t meann = (mean*((Float_t)entries)+value)/((Float_t)entriesn);
301 ((TArrayF *)fPHMean[det])->AddAt(meann,group*fTimeMax+bin);
302 Float_t squaren = ((square*((Float_t)entries))+(value*value))/((Float_t)entriesn);
303 ((TArrayF *)fPHSquares[det])->AddAt(squaren,group*fTimeMax+bin);
311 //__________________________________________________________________________________
312 Bool_t AliTRDCalibraVector::Add(AliTRDCalibraVector *calvect)
315 // Add a other AliTRCalibraVector to this one
318 // Check compatibility
319 if(fNumberBinCharge != calvect->GetNumberBinCharge()) return kFALSE;
320 if(fNumberBinPRF != calvect->GetNumberBinPRF()) return kFALSE;
321 if(fPRFRange != calvect->GetPRFRange()) return kFALSE;
322 if(fTimeMax != calvect->GetTimeMax()) return kFALSE;
323 for(Int_t k = 0; k < 3; k++){
324 if(fDetCha0[k] != calvect->GetDetCha0(k)) return kFALSE;
325 if(fDetCha2[k] != calvect->GetDetCha2(k)) return kFALSE;
328 //printf("pass0!\n");
331 for (Int_t idet = 0; idet < 540; idet++){
333 const TArrayI *phEntriesvect = calvect->GetPHEntries(idet);
334 const TArrayF *phMeanvect = calvect->GetPHMean(idet);
335 const TArrayF *phSquaresvect = calvect->GetPHSquares(idet);
337 const TArrayI *prfEntriesvect = calvect->GetPRFEntries(idet);
338 const TArrayF *prfMeanvect = calvect->GetPRFMean(idet);
339 const TArrayF *prfSquaresvect = calvect->GetPRFSquares(idet);
341 const TArrayI *chEntriesvect = calvect->GetCHEntries(idet);
343 //printf("idet %d!\n",idet);
345 //printf("phEntriesvect %d\n",(Bool_t) phEntriesvect);
346 //printf("phMeanvect %d\n",(Bool_t) phMeanvect);
347 //printf("phSquaresvect %d\n",(Bool_t) phSquaresvect);
349 //printf("prfEntriesvect %d\n",(Bool_t) prfEntriesvect);
350 //printf("prfMeanvect %d\n",(Bool_t) prfMeanvect);
351 //printf("prfSquaresvect %d\n",(Bool_t) prfSquaresvect);
353 //printf("chEntriesvect %d\n",(Bool_t) chEntriesvect);
355 if ( phEntriesvect != 0x0 ){
357 fPHEntries[idet] = ((TArrayI *)GetPHEntries(idet,kTRUE));
358 fPHMean[idet] = ((TArrayF *)GetPHMean(idet,kTRUE));
359 fPHSquares[idet] = ((TArrayF *)GetPHSquares(idet,kTRUE));
360 Int_t total = ((TArrayI *)fPHEntries[idet])->GetSize();
362 for(Int_t k = 0; k < total; k++){
363 Int_t entries = ((TArrayI *)fPHEntries[idet])->At(k);
364 Float_t mean = ((TArrayF *)fPHMean[idet])->At(k);
365 Float_t square = ((TArrayF *)fPHSquares[idet])->At(k);
367 Int_t entriesn = entries+((TArrayI *)phEntriesvect)->At(k);
368 if(entriesn <= 0) continue;
369 ((TArrayI *)fPHEntries[idet])->AddAt(entriesn,k);
370 Float_t meann = (mean*((Float_t)entries)+((TArrayF *)phMeanvect)->At(k)*((Float_t)((TArrayI *)phEntriesvect)->At(k)))/((Float_t)entriesn);
371 ((TArrayF *)fPHMean[idet])->AddAt(meann,k);
372 Float_t sq = ((TArrayF *)phSquaresvect)->At(k)*((Float_t)((TArrayI *)phEntriesvect)->At(k));
373 Float_t squaren = ((square*((Float_t)entries))+sq)/((Float_t)entriesn);
374 ((TArrayF *)fPHSquares[idet])->AddAt(squaren,k);
375 //printf("test ph!\n");
379 if ( prfEntriesvect != 0x0 ){
381 fPRFEntries[idet] = ((TArrayI *)GetPRFEntries(idet,kTRUE));
382 fPRFMean[idet] = ((TArrayF *)GetPRFMean(idet,kTRUE));
383 fPRFSquares[idet] = ((TArrayF *)GetPRFSquares(idet,kTRUE));
384 Int_t total = fPRFEntries[idet]->GetSize();
386 for(Int_t k = 0; k < total; k++){
387 Int_t entries = ((TArrayI *)fPRFEntries[idet])->At(k);
388 Float_t mean = ((TArrayF *)fPRFMean[idet])->At(k);
389 Float_t square = ((TArrayF *)fPRFSquares[idet])->At(k);
391 //printf("entries0 %d\n",entries);
392 //printf("mean0 %f\n",mean);
393 //printf("square0 %f\n",square);
395 //printf("entries1 %d\n",prfEntriesvect->At(k));
396 //printf("mean1 %f\n",prfMeanvect->At(k));
397 //printf("square1 %f\n",prfSquaresvect->At(k));
400 //printf("entries0 size %d\n",fPRFEntries->GetSize());
401 //printf("mean0 size %d\n",fPRFMean->GetSize());
402 //printf("squares0 size %d\n",fPRFSquares->GetSize());
404 //printf("entries1 size %d\n",prfEntriesvect->GetSize());
405 //printf("mean1 size %d\n",prfMeanvect->GetSize());
406 //printf("squares1 size %d\n",prfSquaresvect->GetSize());
409 Int_t entriesn = entries+((TArrayI *)prfEntriesvect)->At(k);
410 if(entriesn <= 0) continue;
411 ((TArrayI *)fPRFEntries[idet])->AddAt(entriesn,k);
412 Float_t meann = (mean*((Float_t)entries)+((TArrayF *)prfMeanvect)->At(k)*((Float_t)((TArrayI *)prfEntriesvect)->At(k)))/((Float_t)entriesn);
413 ((TArrayF *)fPRFMean[idet])->AddAt(meann,k);
414 Float_t sq = ((TArrayF *)prfSquaresvect)->At(k)*((Float_t)((TArrayI *)prfEntriesvect)->At(k));
415 Float_t squaren = ((square*((Float_t)entries))+sq)/((Float_t)entriesn);
416 ((TArrayF *)fPRFSquares[idet])->AddAt(squaren,k);
417 //printf("test prf!\n");
421 if ( chEntriesvect != 0x0 ){
423 fCHEntries[idet] = ((TArrayI *)GetCHEntries(idet,kTRUE));
424 Int_t total = fCHEntries[idet]->GetSize();
425 //if(idet == 180) printf("total %d\n",total);
427 for(Int_t k = 0; k < total; k++){
428 Int_t entries = ((TArrayI *)fCHEntries[idet])->At(k);
429 Int_t entriesn = entries+((TArrayI *)chEntriesvect)->At(k);
430 //if((idet == 180) && ((entries != 0) || (entriesn != 0))) printf("for k %d we have entries %d and entriesn %d\n",k,entries,entriesn);
431 if(entriesn <= 0) continue;
432 ((TArrayI *)fCHEntries[idet])->AddAt(entriesn,k);
434 //printf("test ch!\n");
440 //_____________________________________________________________________________
441 TGraphErrors *AliTRDCalibraVector::ConvertVectorPHTGraphErrors(Int_t det, Int_t group
442 , const Char_t * name)
445 // Convert the fVectorPHMean, fVectorPHSquares and fVectorPHEntries in TGraphErrors
449 fPHEntries[det] = ((TArrayI *)GetPHEntries(det,kTRUE));
450 fPHMean[det] = ((TArrayF *)GetPHMean(det,kTRUE));
451 fPHSquares[det] = ((TArrayF *)GetPHSquares(det,kTRUE));
457 AliTRDCommonParam *parCom = AliTRDCommonParam::Instance();
459 AliInfo("Could not get CommonParam, take the default 10MHz");
461 sf = parCom->GetSamplingFrequency();
469 x = new Double_t[fTimeMax]; // Xaxis
470 y = new Double_t[fTimeMax]; // Sum/entries
471 ex = new Double_t[fTimeMax]; // Nentries
472 ey = new Double_t[fTimeMax]; // Sum of square/nentries
475 Int_t offset = group*fTimeMax;
478 for (Int_t k = 0; k < fTimeMax; k++) {
483 Int_t bin = offset+k;
484 // Fill only if there is more than 0 something
485 if (fPHEntries[det]->At(bin) > 0) {
486 ex[k] = ((TArrayI *)fPHEntries[det])->At(bin);
487 y[k] = ((TArrayF *)fPHMean[det])->At(bin);
488 ey[k] = ((TArrayF *)fPHSquares[det])->At(bin);
492 // Define the TGraphErrors
493 histo = new TGraphErrors(fTimeMax,x,y,ex,ey);
494 histo->SetTitle(name);
498 //_____________________________________________________________________________
499 TGraphErrors *AliTRDCalibraVector::ConvertVectorPRFTGraphErrors(Int_t det, Int_t group
500 , const Char_t * name)
503 // Convert the fVectorPRFMean, fVectorPRFSquares and fVectorPRFEntries in TGraphErrors
507 fPRFEntries[det] = ((TArrayI *)GetPRFEntries(det,kTRUE));
508 fPRFMean[det] = ((TArrayF *)GetPRFMean(det,kTRUE));
509 fPRFSquares[det] = ((TArrayF *)GetPRFSquares(det,kTRUE));
521 x = new Double_t[fNumberBinPRF]; // Xaxis
522 y = new Double_t[fNumberBinPRF]; // Sum/entries
523 ex = new Double_t[fNumberBinPRF]; // Nentries
524 ey = new Double_t[fNumberBinPRF]; // Sum of square/nentries
525 step = (2*TMath::Abs(fPRFRange)) / fNumberBinPRF;
526 min = -TMath::Abs(fPRFRange) + step / 2.0;
527 Int_t offset = group*fNumberBinPRF;
528 //printf("number of total: %d\n",fNumberBinPRF);
530 for (Int_t k = 0; k < fNumberBinPRF; k++) {
535 Int_t bin = offset+k;
536 // Fill only if there is more than 0 something
537 if (fPRFEntries[det]->At(bin) > 0) {
538 ex[k] = ((TArrayF *)fPRFEntries[det])->At(bin);
539 y[k] = ((TArrayF *)fPRFMean[det])->At(bin);
540 ey[k] = ((TArrayF *)fPRFSquares[det])->At(bin);
542 //printf("Number of entries %f for %d\n",ex[k],k);
545 // Define the TGraphErrors
546 histo = new TGraphErrors(fNumberBinPRF,x,y,ex,ey);
547 histo->SetTitle(name);
551 //_____________________________________________________________________________
552 TH1F *AliTRDCalibraVector::ConvertVectorCHHisto(Int_t det, Int_t group
553 , const Char_t * name)
556 // Convert the fVectorCHEntries in TH1F
560 fCHEntries[det] = ((TArrayI *)GetCHEntries(det,kTRUE));
563 TH1F *histo = new TH1F(name,name,fNumberBinCharge,0,300);
565 Int_t offset = group*fNumberBinCharge;
567 for (Int_t k = 0; k < fNumberBinCharge; k++) {
568 Int_t bin = offset+k;
569 histo->SetBinContent(k+1,((TArrayI *)fCHEntries[det])->At(bin));
570 histo->SetBinError(k+1,TMath::Sqrt(TMath::Abs(((TArrayI *)fCHEntries[det])->At(bin))));
576 //_____________________________________________________________________
577 TArrayI* AliTRDCalibraVector::GetPHEntries(Int_t det
578 , Bool_t force) /*FOLD00*/
581 // return pointer to Carge ROC Calibration
582 // if force is true create a new histogram if it doesn't exist allready
584 TArrayI**arr = &fPHEntries[0];
585 return GetEntriesPH(det, arr, force);
587 //_____________________________________________________________________
588 TArrayI* AliTRDCalibraVector::GetPRFEntries(Int_t det
589 , Bool_t force) /*FOLD00*/
592 // return pointer to Carge ROC Calibration
593 // if force is true create a new histogram if it doesn't exist allready
595 TArrayI**arr = &fPRFEntries[0];
596 return GetEntriesPRF(det, arr, force);
598 //_____________________________________________________________________
599 TArrayI* AliTRDCalibraVector::GetCHEntries(Int_t det
600 , Bool_t force) /*FOLD00*/
603 // return pointer to Carge ROC Calibration
604 // if force is true create a new histogram if it doesn't exist allready
606 TArrayI**arr = &fCHEntries[0];
607 return GetEntriesCH(det, arr, force);
609 //_____________________________________________________________________
610 TArrayF* AliTRDCalibraVector::GetPHMean(Int_t det
611 , Bool_t force) /*FOLD00*/
614 // return pointer to Carge ROC Calibration
615 // if force is true create a new histogram if it doesn't exist allready
617 TArrayF**arr = &fPHMean[0];
618 return GetMeanSquaresPH(det, arr, force);
620 //_____________________________________________________________________
621 TArrayF* AliTRDCalibraVector::GetPHSquares(Int_t det
622 , Bool_t force) /*FOLD00*/
625 // return pointer to Carge ROC Calibration
626 // if force is true create a new histogram if it doesn't exist allready
628 TArrayF**arr = &fPHSquares[0];
629 return GetMeanSquaresPH(det, arr, force);
631 //_____________________________________________________________________
632 TArrayF* AliTRDCalibraVector::GetPRFMean(Int_t det
633 , Bool_t force) /*FOLD00*/
636 // return pointer to Carge ROC Calibration
637 // if force is true create a new histogram if it doesn't exist allready
639 TArrayF**arr = &fPRFMean[0];
640 return GetMeanSquaresPRF(det, arr, force);
642 //_____________________________________________________________________
643 TArrayF* AliTRDCalibraVector::GetPRFSquares(Int_t det
644 , Bool_t force) /*FOLD00*/
647 // return pointer to Carge ROC Calibration
648 // if force is true create a new histogram if it doesn't exist allready
650 TArrayF**arr = &fPRFSquares[0];
651 return GetMeanSquaresPRF(det, arr, force);
653 //_____________________________________________________________________
654 TArrayI* AliTRDCalibraVector::GetEntriesCH(Int_t det
656 , Bool_t force) /*FOLD00*/
659 // return pointer to TArrayI Entries
660 // if force is true create a new TArrayI if it doesn't exist allready
662 if ( !force || (((TArrayI *)arr[det])->GetSize()>0))
663 return (TArrayI*)arr[det];
665 // if we are forced and TArrayI doesn't yes exist create it
666 Int_t stack = GetStack(det);
668 if(stack == 2) ngroup = fDetCha2[0]*fNumberBinCharge;
669 else ngroup = fDetCha0[0]*fNumberBinCharge;
671 ((TArrayI *)arr[det])->Set(ngroup);
672 for(Int_t k = 0; k < ngroup; k++){
673 ((TArrayI *)arr[det])->AddAt(0,k);
675 return (TArrayI*)arr[det];
677 //_____________________________________________________________________
678 TArrayI* AliTRDCalibraVector::GetEntriesPRF(Int_t det
680 , Bool_t force) /*FOLD00*/
683 // return pointer to TArrayI Entries
684 // if force is true create a new TArrayI if it doesn't exist allready
686 if ( !force || (((TArrayI *)arr[det])->GetSize()>0))
687 return (TArrayI*)arr[det];
689 // if we are forced and TArrayI doesn't yes exist create it
690 Int_t stack = GetStack(det);
692 if(stack == 2) ngroup = fDetCha2[2]*fNumberBinPRF;
693 else ngroup = fDetCha0[2]*fNumberBinPRF;
695 ((TArrayI *)arr[det])->Set(ngroup);
696 for(Int_t k = 0; k < ngroup; k++){
697 ((TArrayI *)arr[det])->AddAt(0,k);
699 return (TArrayI*)arr[det];
702 //_____________________________________________________________________
703 TArrayI* AliTRDCalibraVector::GetEntriesPH(Int_t det
705 , Bool_t force) /*FOLD00*/
708 // return pointer to TArrayI Entries
709 // if force is true create a new TArrayI if it doesn't exist allready
711 if ( !force || (((TArrayI *)arr[det])->GetSize()>0))
712 return (TArrayI*)arr[det];
714 // if we are forced and TArrayI doesn't yes exist create it
715 Int_t stack = GetStack(det);
717 if(stack == 2) ngroup = fDetCha2[1]*fTimeMax;
718 else ngroup = fDetCha0[1]*fTimeMax;
720 ((TArrayI *)arr[det])->Set(ngroup);
721 for(Int_t k = 0; k < ngroup; k++){
722 ((TArrayI *)arr[det])->AddAt(0,k);
724 return (TArrayI*)arr[det];
727 //_____________________________________________________________________
728 TArrayF* AliTRDCalibraVector::GetMeanSquaresPH(Int_t det
730 , Bool_t force) /*FOLD00*/
733 // return pointer to TArrayF Mean or Squares
734 // if force is true create a new TArrayF if it doesn't exist allready
736 if ( !force || (((TArrayF *)arr[det])->GetSize()>0))
737 return (TArrayF*)arr[det];
739 // if we are forced and TArrayF doesn't yes exist create it
740 Int_t stack = GetStack(det);
742 if(stack == 2) ngroup = fDetCha2[1]*fTimeMax;
743 else ngroup = fDetCha0[1]*fTimeMax;
745 ((TArrayF *)arr[det])->Set(ngroup);
746 for(Int_t k = 0; k < ngroup; k++){
747 ((TArrayF *)arr[det])->AddAt(0.0,k);
749 return ((TArrayF *)arr[det]);
751 //_____________________________________________________________________
752 TArrayF* AliTRDCalibraVector::GetMeanSquaresPRF(Int_t det
754 , Bool_t force) /*FOLD00*/
757 // return pointer to TArrayF Mean or Squares
758 // if force is true create a new TArrayF if it doesn't exist allready
760 if ( !force || (((TArrayF *)arr[det])->GetSize()>0))
763 // if we are forced and TArrayF doesn't yes exist create it
764 Int_t stack = GetStack(det);
766 if(stack == 2) ngroup = fDetCha2[2]*fNumberBinPRF;
767 else ngroup = fDetCha0[2]*fNumberBinPRF;
769 ((TArrayF *)arr[det])->Set(ngroup);
770 for(Int_t k = 0; k < ngroup; k++){
771 ((TArrayF *)arr[det])->AddAt(0.0,k);
773 return ((TArrayF *)arr[det]);
775 //_____________________________________________________________________________
776 Int_t AliTRDCalibraVector::GetStack(Int_t d) const
779 // Reconstruct the stack number from the detector number
782 return ((Int_t) (d % 30) / 6);