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"
41 #include "AliTRDarrayF.h"
42 #include "AliTRDarrayI.h"
44 ClassImp(AliTRDCalibraVector)
46 //______________________________________________________________________________________
47 AliTRDCalibraVector::AliTRDCalibraVector()
60 ,fVectorPHSquares(540)
61 ,fVectorPHEntries(540)
62 ,fVectorCHEntries(540)
64 ,fVectorPRFSquares(540)
65 ,fVectorPRFEntries(540)
72 // Default constructor
75 for(Int_t k = 0; k < 3; k++){
81 //______________________________________________________________________________________
82 AliTRDCalibraVector::AliTRDCalibraVector(const AliTRDCalibraVector &c)
95 ,fVectorPHSquares(540)
96 ,fVectorPHEntries(540)
97 ,fVectorCHEntries(540)
99 ,fVectorPRFSquares(540)
100 ,fVectorPRFEntries(540)
101 ,fNumberBinCharge(c.fNumberBinCharge)
102 ,fNumberBinPRF(c.fNumberBinPRF)
103 ,fTimeMax(c.fTimeMax)
104 ,fPRFRange(c.fPRFRange)
109 for (Int_t idet = 0; idet < 540; idet++){
111 const AliTRDarrayI *phEntries = (AliTRDarrayI*)c.fVectorPHEntries.UncheckedAt(idet);
112 const AliTRDarrayF *phMean = (AliTRDarrayF*)c.fVectorPHMean.UncheckedAt(idet);
113 const AliTRDarrayF *phSquares = (AliTRDarrayF*)c.fVectorPHSquares.UncheckedAt(idet);
115 const AliTRDarrayI *prfEntries = (AliTRDarrayI*)c.fVectorPRFEntries.UncheckedAt(idet);
116 const AliTRDarrayF *prfMean = (AliTRDarrayF*)c.fVectorPRFMean.UncheckedAt(idet);
117 const AliTRDarrayF *prfSquares = (AliTRDarrayF*)c.fVectorPRFSquares.UncheckedAt(idet);
119 const AliTRDarrayI *chEntries = (AliTRDarrayI*)c.fVectorCHEntries.UncheckedAt(idet);
121 if ( phEntries != 0x0 ) fVectorPHEntries.AddAt(new AliTRDarrayI(*phEntries), idet);
122 if ( phMean != 0x0 ) fVectorPHMean.AddAt(new AliTRDarrayF(*phMean), idet);
123 if ( phSquares != 0x0 ) fVectorPHSquares.AddAt(new AliTRDarrayF(*phSquares), idet);
125 if ( prfEntries != 0x0 ) fVectorPRFEntries.AddAt(new AliTRDarrayI(*prfEntries), idet);
126 if ( prfMean != 0x0 ) fVectorPRFMean.AddAt(new AliTRDarrayF(*prfMean), idet);
127 if ( prfSquares != 0x0 ) fVectorPRFSquares.AddAt(new AliTRDarrayF(*prfSquares), idet);
130 if ( chEntries != 0x0 ) fVectorCHEntries.AddAt(new AliTRDarrayI(*chEntries), idet);
134 for(Int_t k = 0; k < 3; k++){
135 fDetCha0[k] = c.fDetCha0[k];
136 fDetCha2[k] = c.fDetCha2[k];
139 fVectorPHEntries.SetName(c.fVectorPHEntries.GetName());
140 fVectorCHEntries.SetName(c.fVectorCHEntries.GetName());
141 fVectorPRFEntries.SetName(c.fVectorPRFEntries.GetName());
144 //_____________________________________________________________________
145 AliTRDCalibraVector& AliTRDCalibraVector::operator = (const AliTRDCalibraVector &source)
148 // assignment operator
150 if (&source == this) return *this;
151 new (this) AliTRDCalibraVector(source);
155 //____________________________________________________________________________________
156 AliTRDCalibraVector::~AliTRDCalibraVector()
159 // AliTRDCalibraVector destructor
163 //_____________________________________________________________________________
164 Int_t AliTRDCalibraVector::SearchBin(Float_t value, Int_t i) const
172 Float_t fbinmax = value;
173 Int_t fNumberOfBin = -1;
180 fNumberOfBin = fNumberBinCharge;
184 fbinmax = TMath::Abs(fPRFRange);
185 fbinmin = -TMath::Abs(fPRFRange);
186 fNumberOfBin = fNumberBinPRF;
194 if ((value >= fbinmax) ||
199 reponse = (Int_t) ((fNumberOfBin*(value-fbinmin)) / (fbinmax-fbinmin));
205 //_____________________________________________________________________________
206 Bool_t AliTRDCalibraVector::UpdateVectorCH(Int_t det, Int_t group, Float_t value)
209 // Fill the vector CH
213 Int_t bin = SearchBin(value,0);
221 if(fDetectorCH != det){
222 fCHEntries = ((AliTRDarrayI *)GetCHEntries(det,kTRUE));
225 Int_t entries = fCHEntries->At(group*fNumberBinCharge+bin);
227 Int_t entriesn = entries+1;
228 fCHEntries->AddAt(entriesn,group*fNumberBinCharge+bin);
236 //_____________________________________________________________________________
237 Bool_t AliTRDCalibraVector::UpdateVectorPRF(Int_t det, Int_t group, Float_t x, Float_t y)
240 // Fill the vector PRF
244 Int_t bin = SearchBin(x,2);
251 if(fDetectorPRF != det){
252 fPRFEntries = ((AliTRDarrayI *)GetPRFEntries(det,kTRUE));
253 fPRFMean = ((AliTRDarrayF *)GetPRFMean(det,kTRUE));
254 fPRFSquares = ((AliTRDarrayF *)GetPRFSquares(det,kTRUE));
257 Int_t entries = fPRFEntries->At(group*fNumberBinPRF+bin);
258 Float_t mean = fPRFMean->At(group*fNumberBinPRF+bin);
259 Float_t square = fPRFSquares->At(group*fNumberBinPRF+bin);
261 Int_t entriesn = entries+1;
262 fPRFEntries->AddAt(entriesn,group*fNumberBinPRF+bin);
263 Float_t meann = (mean*((Float_t)entries)+y)/((Float_t)entriesn);
264 fPRFMean->AddAt(meann,group*fNumberBinPRF+bin);
265 Float_t squaren = ((square*((Float_t)entries))+(y*y))/((Float_t)entriesn);
266 fPRFSquares->AddAt(squaren,group*fNumberBinPRF+bin);
274 //_____________________________________________________________________________
275 Bool_t AliTRDCalibraVector::UpdateVectorPH(Int_t det, Int_t group, Int_t time, Float_t value)
278 // Fill the vector PH
290 if(fDetectorPH != det){
291 fPHEntries = ((AliTRDarrayI *)GetPHEntries(det,kTRUE));
292 fPHMean = ((AliTRDarrayF *)GetPHMean(det,kTRUE));
293 fPHSquares = ((AliTRDarrayF *)GetPHSquares(det,kTRUE));
296 Int_t entries = fPHEntries->At(group*fTimeMax+bin);
297 Float_t mean = fPHMean->At(group*fTimeMax+bin);
298 Float_t square = fPHSquares->At(group*fTimeMax+bin);
300 Int_t entriesn = entries+1;
301 fPHEntries->AddAt(entriesn,group*fTimeMax+bin);
302 Float_t meann = (mean*((Float_t)entries)+value)/((Float_t)entriesn);
303 fPHMean->AddAt(meann,group*fTimeMax+bin);
304 Float_t squaren = ((square*((Float_t)entries))+(value*value))/((Float_t)entriesn);
305 fPHSquares->AddAt(squaren,group*fTimeMax+bin);
313 //__________________________________________________________________________________
314 Bool_t AliTRDCalibraVector::Add(AliTRDCalibraVector *calvect)
317 // Add a other AliTRCalibraVector to this one
320 // Check compatibility
321 if(fNumberBinCharge != calvect->GetNumberBinCharge()) return kFALSE;
322 if(fNumberBinPRF != calvect->GetNumberBinPRF()) return kFALSE;
323 if(fPRFRange != calvect->GetPRFRange()) return kFALSE;
324 if(fTimeMax != calvect->GetTimeMax()) return kFALSE;
325 for(Int_t k = 0; k < 3; k++){
326 if(fDetCha0[k] != calvect->GetDetCha0(k)) return kFALSE;
327 if(fDetCha2[k] != calvect->GetDetCha2(k)) return kFALSE;
330 //printf("pass0!\n");
333 for (Int_t idet = 0; idet < 540; idet++){
335 const AliTRDarrayI *phEntriesvect = calvect->GetPHEntries(idet);
336 const AliTRDarrayF *phMeanvect = calvect->GetPHMean(idet);
337 const AliTRDarrayF *phSquaresvect = calvect->GetPHSquares(idet);
339 const AliTRDarrayI *prfEntriesvect = calvect->GetPRFEntries(idet);
340 const AliTRDarrayF *prfMeanvect = calvect->GetPRFMean(idet);
341 const AliTRDarrayF *prfSquaresvect = calvect->GetPRFSquares(idet);
343 const AliTRDarrayI *chEntriesvect = calvect->GetCHEntries(idet);
345 //printf("idet %d!\n",idet);
347 //printf("phEntriesvect %d\n",(Bool_t) phEntriesvect);
348 //printf("phMeanvect %d\n",(Bool_t) phMeanvect);
349 //printf("phSquaresvect %d\n",(Bool_t) phSquaresvect);
351 //printf("prfEntriesvect %d\n",(Bool_t) prfEntriesvect);
352 //printf("prfMeanvect %d\n",(Bool_t) prfMeanvect);
353 //printf("prfSquaresvect %d\n",(Bool_t) prfSquaresvect);
355 //printf("chEntriesvect %d\n",(Bool_t) chEntriesvect);
357 if ( phEntriesvect != 0x0 ){
359 fPHEntries = ((AliTRDarrayI *)GetPHEntries(idet,kTRUE));
360 fPHMean = ((AliTRDarrayF *)GetPHMean(idet,kTRUE));
361 fPHSquares = ((AliTRDarrayF *)GetPHSquares(idet,kTRUE));
362 Int_t total = fPHEntries->GetSize();
364 for(Int_t k = 0; k < total; k++){
365 Int_t entries = fPHEntries->At(k);
366 Float_t mean = fPHMean->At(k);
367 Float_t square = fPHSquares->At(k);
369 Int_t entriesn = entries+phEntriesvect->At(k);
370 if(entriesn <= 0) continue;
371 fPHEntries->AddAt(entriesn,k);
372 Float_t meann = (mean*((Float_t)entries)+phMeanvect->At(k)*((Float_t)phEntriesvect->At(k)))/((Float_t)entriesn);
373 fPHMean->AddAt(meann,k);
374 Float_t sq = phSquaresvect->At(k)*((Float_t)phEntriesvect->At(k));
375 Float_t squaren = ((square*((Float_t)entries))+sq)/((Float_t)entriesn);
376 fPHSquares->AddAt(squaren,k);
377 //printf("test ph!\n");
381 if ( prfEntriesvect != 0x0 ){
383 fPRFEntries = ((AliTRDarrayI *)GetPRFEntries(idet,kTRUE));
384 fPRFMean = ((AliTRDarrayF *)GetPRFMean(idet,kTRUE));
385 fPRFSquares = ((AliTRDarrayF *)GetPRFSquares(idet,kTRUE));
386 Int_t total = fPRFEntries->GetSize();
388 for(Int_t k = 0; k < total; k++){
389 Int_t entries = fPRFEntries->At(k);
390 Float_t mean = fPRFMean->At(k);
391 Float_t square = fPRFSquares->At(k);
393 //printf("entries0 %d\n",entries);
394 //printf("mean0 %f\n",mean);
395 //printf("square0 %f\n",square);
397 //printf("entries1 %d\n",prfEntriesvect->At(k));
398 //printf("mean1 %f\n",prfMeanvect->At(k));
399 //printf("square1 %f\n",prfSquaresvect->At(k));
402 //printf("entries0 size %d\n",fPRFEntries->GetSize());
403 //printf("mean0 size %d\n",fPRFMean->GetSize());
404 //printf("squares0 size %d\n",fPRFSquares->GetSize());
406 //printf("entries1 size %d\n",prfEntriesvect->GetSize());
407 //printf("mean1 size %d\n",prfMeanvect->GetSize());
408 //printf("squares1 size %d\n",prfSquaresvect->GetSize());
411 Int_t entriesn = entries+prfEntriesvect->At(k);
412 if(entriesn <= 0) continue;
413 fPRFEntries->AddAt(entriesn,k);
414 Float_t meann = (mean*((Float_t)entries)+prfMeanvect->At(k)*((Float_t)prfEntriesvect->At(k)))/((Float_t)entriesn);
415 fPRFMean->AddAt(meann,k);
416 Float_t sq = prfSquaresvect->At(k)*((Float_t)prfEntriesvect->At(k));
417 Float_t squaren = ((square*((Float_t)entries))+sq)/((Float_t)entriesn);
418 fPRFSquares->AddAt(squaren,k);
419 //printf("test prf!\n");
423 if ( chEntriesvect != 0x0 ){
425 fCHEntries = ((AliTRDarrayI *)GetCHEntries(idet,kTRUE));
426 Int_t total = fCHEntries->GetSize();
427 //if(idet == 180) printf("total %d\n",total);
429 for(Int_t k = 0; k < total; k++){
430 Int_t entries = fCHEntries->At(k);
431 Int_t entriesn = entries+chEntriesvect->At(k);
432 //if((idet == 180) && ((entries != 0) || (entriesn != 0))) printf("for k %d we have entries %d and entriesn %d\n",k,entries,entriesn);
433 if(entriesn <= 0) continue;
434 fCHEntries->AddAt(entriesn,k);
436 //printf("test ch!\n");
442 //_____________________________________________________________________________
443 TGraphErrors *AliTRDCalibraVector::ConvertVectorPHTGraphErrors(Int_t det, Int_t group
444 , const Char_t * name)
447 // Convert the fVectorPHMean, fVectorPHSquares and fVectorPHEntries in TGraphErrors
451 fPHEntries = ((AliTRDarrayI *)GetPHEntries(det,kTRUE));
452 fPHMean = ((AliTRDarrayF *)GetPHMean(det,kTRUE));
453 fPHSquares = ((AliTRDarrayF *)GetPHSquares(det,kTRUE));
459 AliTRDCommonParam *parCom = AliTRDCommonParam::Instance();
461 AliInfo("Could not get CommonParam, take the default 10MHz");
463 sf = parCom->GetSamplingFrequency();
471 x = new Double_t[fTimeMax]; // Xaxis
472 y = new Double_t[fTimeMax]; // Sum/entries
473 ex = new Double_t[fTimeMax]; // Nentries
474 ey = new Double_t[fTimeMax]; // Sum of square/nentries
477 Int_t offset = group*fTimeMax;
480 for (Int_t k = 0; k < fTimeMax; k++) {
485 Int_t bin = offset+k;
486 // Fill only if there is more than 0 something
487 if (fPHEntries->At(bin) > 0) {
488 ex[k] = fPHEntries->At(bin);
489 y[k] = fPHMean->At(bin);
490 ey[k] = fPHSquares->At(bin);
494 // Define the TGraphErrors
495 histo = new TGraphErrors(fTimeMax,x,y,ex,ey);
496 histo->SetTitle(name);
500 //_____________________________________________________________________________
501 TGraphErrors *AliTRDCalibraVector::ConvertVectorPRFTGraphErrors(Int_t det, Int_t group
502 , const Char_t * name)
505 // Convert the fVectorPRFMean, fVectorPRFSquares and fVectorPRFEntries in TGraphErrors
509 fPRFEntries = ((AliTRDarrayI *)GetPRFEntries(det,kTRUE));
510 fPRFMean = ((AliTRDarrayF *)GetPRFMean(det,kTRUE));
511 fPRFSquares = ((AliTRDarrayF *)GetPRFSquares(det,kTRUE));
523 x = new Double_t[fNumberBinPRF]; // Xaxis
524 y = new Double_t[fNumberBinPRF]; // Sum/entries
525 ex = new Double_t[fNumberBinPRF]; // Nentries
526 ey = new Double_t[fNumberBinPRF]; // Sum of square/nentries
527 step = (2*TMath::Abs(fPRFRange)) / fNumberBinPRF;
528 min = -TMath::Abs(fPRFRange) + step / 2.0;
529 Int_t offset = group*fNumberBinPRF;
530 //printf("number of total: %d\n",fNumberBinPRF);
532 for (Int_t k = 0; k < fNumberBinPRF; k++) {
537 Int_t bin = offset+k;
538 // Fill only if there is more than 0 something
539 if (fPRFEntries->At(bin) > 0) {
540 ex[k] = fPRFEntries->At(bin);
541 y[k] = fPRFMean->At(bin);
542 ey[k] = fPRFSquares->At(bin);
544 //printf("Number of entries %f for %d\n",ex[k],k);
547 // Define the TGraphErrors
548 histo = new TGraphErrors(fNumberBinPRF,x,y,ex,ey);
549 histo->SetTitle(name);
553 //_____________________________________________________________________________
554 TH1F *AliTRDCalibraVector::ConvertVectorCHHisto(Int_t det, Int_t group
555 , const Char_t * name)
558 // Convert the fVectorCHEntries in TH1F
562 fCHEntries = ((AliTRDarrayI *)GetCHEntries(det,kTRUE));
565 TH1F *histo = new TH1F(name,name,fNumberBinCharge,0,300);
567 Int_t offset = group*fNumberBinCharge;
569 for (Int_t k = 0; k < fNumberBinCharge; k++) {
570 Int_t bin = offset+k;
571 histo->SetBinContent(k+1,fCHEntries->At(bin));
572 histo->SetBinError(k+1,TMath::Sqrt(TMath::Abs(fCHEntries->At(bin))));
578 //_____________________________________________________________________
579 AliTRDarrayI* AliTRDCalibraVector::GetPHEntries(Int_t det
580 , Bool_t force) /*FOLD00*/
583 // return pointer to Carge ROC Calibration
584 // if force is true create a new histogram if it doesn't exist allready
586 TObjArray *arr = &fVectorPHEntries;
587 return GetEntriesPH(det, arr, force);
589 //_____________________________________________________________________
590 AliTRDarrayI* AliTRDCalibraVector::GetPRFEntries(Int_t det
591 , Bool_t force) /*FOLD00*/
594 // return pointer to Carge ROC Calibration
595 // if force is true create a new histogram if it doesn't exist allready
597 TObjArray *arr = &fVectorPRFEntries;
598 return GetEntriesPRF(det, arr, force);
600 //_____________________________________________________________________
601 AliTRDarrayI* AliTRDCalibraVector::GetCHEntries(Int_t det
602 , Bool_t force) /*FOLD00*/
605 // return pointer to Carge ROC Calibration
606 // if force is true create a new histogram if it doesn't exist allready
608 TObjArray *arr = &fVectorCHEntries;
609 return GetEntriesCH(det, arr, force);
611 //_____________________________________________________________________
612 AliTRDarrayF* AliTRDCalibraVector::GetPHMean(Int_t det
613 , Bool_t force) /*FOLD00*/
616 // return pointer to Carge ROC Calibration
617 // if force is true create a new histogram if it doesn't exist allready
619 TObjArray *arr = &fVectorPHMean;
620 return GetMeanSquaresPH(det, arr, force);
622 //_____________________________________________________________________
623 AliTRDarrayF* AliTRDCalibraVector::GetPHSquares(Int_t det
624 , Bool_t force) /*FOLD00*/
627 // return pointer to Carge ROC Calibration
628 // if force is true create a new histogram if it doesn't exist allready
630 TObjArray *arr = &fVectorPHSquares;
631 return GetMeanSquaresPH(det, arr, force);
633 //_____________________________________________________________________
634 AliTRDarrayF* AliTRDCalibraVector::GetPRFMean(Int_t det
635 , Bool_t force) /*FOLD00*/
638 // return pointer to Carge ROC Calibration
639 // if force is true create a new histogram if it doesn't exist allready
641 TObjArray *arr = &fVectorPRFMean;
642 return GetMeanSquaresPRF(det, arr, force);
644 //_____________________________________________________________________
645 AliTRDarrayF* AliTRDCalibraVector::GetPRFSquares(Int_t det
646 , Bool_t force) /*FOLD00*/
649 // return pointer to Carge ROC Calibration
650 // if force is true create a new histogram if it doesn't exist allready
652 TObjArray *arr = &fVectorPRFSquares;
653 return GetMeanSquaresPRF(det, arr, force);
655 //_____________________________________________________________________
656 AliTRDarrayI* AliTRDCalibraVector::GetEntriesCH(Int_t det
658 , Bool_t force) /*FOLD00*/
661 // return pointer to AliTRDarrayI Entries
662 // if force is true create a new AliTRDarrayI if it doesn't exist allready
664 if ( !force || arr->UncheckedAt(det) )
665 return (AliTRDarrayI*)arr->UncheckedAt(det);
667 // if we are forced and AliTRDarrayI doesn't yes exist create it
668 AliTRDarrayI *croc = new AliTRDarrayI();
669 Int_t chamber = GetChamber(det);
671 if(chamber == 2) ngroup = fDetCha2[0]*fNumberBinCharge;
672 else ngroup = fDetCha0[0]*fNumberBinCharge;
674 croc->Expand(ngroup);
675 for(Int_t k = 0; k < ngroup; k++){
678 arr->AddAt(croc,det);
681 //_____________________________________________________________________
682 AliTRDarrayI* AliTRDCalibraVector::GetEntriesPRF(Int_t det
684 , Bool_t force) /*FOLD00*/
687 // return pointer to AliTRDarrayI Entries
688 // if force is true create a new AliTRDarrayI if it doesn't exist allready
690 if ( !force || arr->UncheckedAt(det) )
691 return (AliTRDarrayI*)arr->UncheckedAt(det);
693 // if we are forced and AliTRDarrayI doesn't yes exist create it
694 AliTRDarrayI *croc = new AliTRDarrayI();
695 Int_t chamber = GetChamber(det);
697 if(chamber == 2) ngroup = fDetCha2[2]*fNumberBinPRF;
698 else ngroup = fDetCha0[2]*fNumberBinPRF;
700 croc->Expand(ngroup);
701 for(Int_t k = 0; k < ngroup; k++){
704 arr->AddAt(croc,det);
708 //_____________________________________________________________________
709 AliTRDarrayI* AliTRDCalibraVector::GetEntriesPH(Int_t det
711 , Bool_t force) /*FOLD00*/
714 // return pointer to AliTRDarrayI Entries
715 // if force is true create a new AliTRDarrayI if it doesn't exist allready
717 if ( !force || arr->UncheckedAt(det) )
718 return (AliTRDarrayI*)arr->UncheckedAt(det);
720 // if we are forced and AliTRDarrayI doesn't yes exist create it
721 AliTRDarrayI *croc = new AliTRDarrayI();
722 Int_t chamber = GetChamber(det);
724 if(chamber == 2) ngroup = fDetCha2[1]*fTimeMax;
725 else ngroup = fDetCha0[1]*fTimeMax;
727 croc->Expand(ngroup);
728 for(Int_t k = 0; k < ngroup; k++){
731 arr->AddAt(croc,det);
735 //_____________________________________________________________________
736 AliTRDarrayF* AliTRDCalibraVector::GetMeanSquaresPH(Int_t det
738 , Bool_t force) /*FOLD00*/
741 // return pointer to AliTRDarrayF Mean or Squares
742 // if force is true create a new AliTRDarrayF if it doesn't exist allready
744 if ( !force || arr->UncheckedAt(det) )
745 return (AliTRDarrayF*)arr->UncheckedAt(det);
747 // if we are forced and AliTRDarrayF doesn't yes exist create it
748 AliTRDarrayF *croc = new AliTRDarrayF();
749 Int_t chamber = GetChamber(det);
751 if(chamber == 2) ngroup = fDetCha2[1]*fTimeMax;
752 else ngroup = fDetCha0[1]*fTimeMax;
754 croc->Expand(ngroup);
755 for(Int_t k = 0; k < ngroup; k++){
758 arr->AddAt(croc,det);
761 //_____________________________________________________________________
762 AliTRDarrayF* AliTRDCalibraVector::GetMeanSquaresPRF(Int_t det
764 , Bool_t force) /*FOLD00*/
767 // return pointer to AliTRDarrayF Mean or Squares
768 // if force is true create a new AliTRDarrayF if it doesn't exist allready
770 if ( !force || arr->UncheckedAt(det) )
771 return (AliTRDarrayF*)arr->UncheckedAt(det);
773 // if we are forced and AliTRDarrayF doesn't yes exist create it
774 AliTRDarrayF *croc = new AliTRDarrayF();
775 Int_t chamber = GetChamber(det);
777 if(chamber == 2) ngroup = fDetCha2[2]*fNumberBinPRF;
778 else ngroup = fDetCha0[2]*fNumberBinPRF;
780 croc->Expand(ngroup);
781 for(Int_t k = 0; k < ngroup; k++){
784 arr->AddAt(croc,det);
787 //_____________________________________________________________________________
788 Int_t AliTRDCalibraVector::GetChamber(Int_t d) const
791 // Reconstruct the chamber number from the detector number
794 return ((Int_t) (d % 30) / 6);