TRD/AliTRDCalibPadStatus.cxx

   1 /**************************************************************************
   2  * Copyright(c) 1998-1999, ALICE Experiment at CERN, All rights reserved. *
   3  *                                                                        *
   4  * Author: The ALICE Off-line Project.                                    *
   5  * Contributors are mentioned in the code where appropriate.              *
   6  *                                                                        *
   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  **************************************************************************/
  15
  16 /* $Id$ */
  17
  18 /*
  19  example: fill pedestal with gausschen noise
  20  AliTRDCalibPadStatus ped;
  21  ped.TestEvent(numberofevent);
  22  // Method without histo
  23  //////////////////////////
  24  ped.Analyse();
  25  //Create the histo of the AliTRDCalROC
  26  TH2F * histo2dm = ped.GetCalRocMean(0,kFALSE)->MakeHisto2D();
  27  histo2dm->Scale(10.0);
  28  TH1F * histo1dm = ped.GetCalRocMean(0,kFALSE)->MakeHisto1D();
  29  histo1dm->Scale(10.0);
  30  TH2F * histo2ds = ped.GetCalRocSquares(0,kFALSE)->MakeHisto2D();
  31  histo2ds->Scale(10.0);
  32  TH1F * histo1ds = ped.GetCalRocSquares(0,kFALSE)->MakeHisto1D();
  33  histo1ds->Scale(10.0)
  34  //Draw output;
  35  TCanvas* c1 = new TCanvas;
  36  c1->Divide(2,2);
  37  c1->cd(1);
  38  histo2dm->Draw("colz");
  39  c1->cd(2);
  40  histo1dm->Draw();
  41  c1->cd(3);
  42  histo2ds->Draw("colz");
  43  c1->cd(4);
  44  histo1ds->Draw();
  45 // Method with histo
  46 /////////////////////////
  47 ped.AnalyseHisto();
  48 //Take the histo
  49 TH1F * histo = ped.GetHisto(31);
  50 histo->SetEntries(1);
  51 histo->Draw();
  52
  53 */
  54
  55 //Root includes
  56 #include <TObjArray.h>
  57 #include <TH1F.h>
  58 #include <TH2F.h>
  59 #include <TString.h>
  60 #include <TMath.h>
  61 #include <TF1.h>
  62 #include <TRandom.h>
  63 #include <TDirectory.h>
  64 #include <TFile.h>
  65 #include <AliMathBase.h>
  66 #include "TTreeStream.h"
  67
  68 //AliRoot includes
  69 #include "AliRawReader.h"
  70 #include "AliRawReaderRoot.h"
  71 #include "AliRawReaderDate.h"
  72
  73 #include "AliTRDRawStream.h"
  74 #include "AliTRDarrayF.h"
  75 #include "AliTRDgeometry.h"
  76 #include "./Cal/AliTRDCalROC.h"
  77 #include "./Cal/AliTRDCalPadStatus.h"
  78 #include "./Cal/AliTRDCalSingleChamberStatus.h"
  79
  80 #ifdef ALI_DATE
  81 #include "event.h"
  82 #endif
  83
  84 //header file
  85 #include "AliTRDCalibPadStatus.h"
  86
  87 ClassImp(AliTRDCalibPadStatus) /*FOLD00*/
  88
  89 //_____________________________________________________________________
  90 AliTRDCalibPadStatus::AliTRDCalibPadStatus() : /*FOLD00*/
  91   TObject(),
  92   fGeo(0),
  93   fAdcMin(0),
  94   fAdcMax(20),
  95   fDetector(-1),
  96   fNumberOfTimeBins(0),
  97   fCalArrayEntries(540),
  98   fCalArrayMean(540),
  99   fCalArraySquares(540),
 100   fCalRocArrayMean(540),
 101   fCalRocArrayRMS(540),
 102   fHistoArray(540),
 103   fCalEntries(0x0),
 104   fCalMean(0x0),
 105   fCalSquares(0x0)
 106 {
 107     //
 108     // default constructor
 109     //
 110
 111   fGeo = new AliTRDgeometry();
 112
 113 }
 114
 115 //_____________________________________________________________________
 116 AliTRDCalibPadStatus::AliTRDCalibPadStatus(const AliTRDCalibPadStatus &ped) : /*FOLD00*/
 117   TObject(ped),
 118   fGeo(0),
 119   fAdcMin(ped.GetAdcMin()),
 120   fAdcMax(ped.GetAdcMax()),
 121   fDetector(ped.fDetector),
 122   fNumberOfTimeBins(ped.fNumberOfTimeBins),
 123   fCalArrayEntries(540),
 124   fCalArrayMean(540),
 125   fCalArraySquares(540),
 126   fCalRocArrayMean(540),
 127   fCalRocArrayRMS(540),
 128   fHistoArray(540),
 129   fCalEntries(0x0),
 130   fCalMean(0x0),
 131   fCalSquares(0x0)
 132 {
 133     //
 134     // copy constructor
 135     //
 136     for (Int_t idet = 0; idet < 540; idet++){
 137         const AliTRDarrayF *calEntries  = (AliTRDarrayF*)ped.fCalArrayEntries.UncheckedAt(idet);
 138         const AliTRDarrayF *calMean     = (AliTRDarrayF*)ped.fCalArrayMean.UncheckedAt(idet);
 139         const AliTRDarrayF *calSquares  = (AliTRDarrayF*)ped.fCalArraySquares.UncheckedAt(idet);
 140         const AliTRDCalROC *calRocMean  = (AliTRDCalROC*)ped.fCalRocArrayMean.UncheckedAt(idet);
 141         const AliTRDCalROC *calRocRMS   = (AliTRDCalROC*)ped.fCalRocArrayRMS.UncheckedAt(idet);
 142         const TH2F         *hped        = (TH2F*)ped.fHistoArray.UncheckedAt(idet);
 143
 144         if ( calEntries != 0x0 ) fCalArrayEntries.AddAt(new AliTRDarrayF(*calEntries), idet);
 145         if ( calMean != 0x0 )    fCalArrayMean.AddAt(new AliTRDarrayF(*calMean), idet);
 146         if ( calSquares != 0x0 ) fCalArraySquares.AddAt(new AliTRDarrayF(*calSquares), idet);
 147         if ( calRocMean != 0x0 ) fCalRocArrayMean.AddAt(new AliTRDCalROC(*calRocMean), idet);
 148         if ( calRocRMS != 0x0 )  fCalRocArrayRMS.AddAt(new AliTRDCalROC(*calRocRMS), idet);
 149
 150         if ( hped != 0x0 ){
 151           TH2F *hNew = new TH2F(*hped);
 152           hNew->SetDirectory(0);
 153           fHistoArray.AddAt(hNew,idet);
 154         }
 155
 156     }
 157
 158     if (fGeo) {
 159       delete fGeo;
 160     }
 161     fGeo = new AliTRDgeometry();
 162
 163 }
 164
 165 //_____________________________________________________________________
 166 AliTRDCalibPadStatus& AliTRDCalibPadStatus::operator = (const  AliTRDCalibPadStatus &source)
 167 {
 168   //
 169   // assignment operator
 170   //
 171   if (&source == this) return *this;
 172   new (this) AliTRDCalibPadStatus(source);
 173
 174   return *this;
 175 }
 176
 177 //_____________________________________________________________________
 178 AliTRDCalibPadStatus::~AliTRDCalibPadStatus() /*FOLD00*/
 179 {
 180   //
 181   // destructor
 182   //
 183
 184   if (fGeo) {
 185     delete fGeo;
 186   }
 187
 188 }
 189
 190 //_____________________________________________________________________
 191 Int_t AliTRDCalibPadStatus::Update(const Int_t icdet, /*FOLD00*/
 192                                    const Int_t icRow,
 193                                    const Int_t icCol,
 194                                    const Int_t csignal,
 195                                    const Int_t crowMax)
 196 {
 197     //
 198     // Signal filling methode
 199     //
 200     if ( (csignal>fAdcMax) || (csignal<fAdcMin)   ) return 0;
 201
 202     if(fDetector != icdet){
 203       fCalEntries = ((AliTRDarrayF *)GetCalEntries(icdet,kTRUE));
 204       fCalMean    = ((AliTRDarrayF *)GetCalMean(icdet,kTRUE));
 205       fCalSquares = ((AliTRDarrayF *)GetCalSquares(icdet,kTRUE));
 206     }
 207
 208     Float_t entries  = fCalEntries->At(icRow+icCol*crowMax);
 209     Float_t mean     = fCalMean->At(icRow+icCol*crowMax);
 210     Float_t squares  = fCalSquares->At(icRow+icCol*crowMax);
 211
 212     Float_t entriesn = entries+1.0;
 213     fCalEntries->AddAt(entriesn,(icRow+icCol*crowMax));
 214     Float_t meann = (mean*entries+((Float_t)(csignal+0.5)))/entriesn;
 215     fCalMean->AddAt(meann,icRow+icCol*crowMax);
 216     Float_t squaresn = ((squares*entries)+(((Float_t)(csignal+0.5))*((Float_t)(csignal+0.5))))/entriesn;
 217     fCalSquares->AddAt(squaresn,icRow+icCol*crowMax);
 218
 219     //printf("icdet %d, icRow %d, icCol %d, csignal %d, crowMax %d\n",icdet,icRow,icCol,csignal,crowMax);
 220     //printf("entries %f, mean %f, squares %f\n",entriesn,meann,squaresn);
 221
 222     fDetector = icdet;
 223
 224     return 0;
 225 }
 226 //_____________________________________________________________________
 227 Int_t AliTRDCalibPadStatus::UpdateHisto(const Int_t icdet, /*FOLD00*/
 228                                        const Int_t icRow,
 229                                        const Int_t icCol,
 230                                        const Int_t csignal,
 231                                        const Int_t crowMax)
 232 {
 233     //
 234     // Signal filling methode
 235     //
 236   Int_t nbchannel = icRow+icCol*crowMax;
 237
 238   // fast filling methode.
 239   // Attention: the entry counter of the histogram is not increased
 240   //            this means that e.g. the colz draw option gives an empty plot
 241   Int_t bin = 0;
 242   if ( !(((Int_t)csignal>fAdcMax ) || ((Int_t)csignal<fAdcMin)) )
 243     bin = (nbchannel+1)*(fAdcMax-fAdcMin+2)+((Int_t)csignal-fAdcMin+1);
 244
 245   GetHisto(icdet,kTRUE)->GetArray()[bin]++;
 246
 247   return 0;
 248 }
 249 //_____________________________________________________________________
 250 Bool_t AliTRDCalibPadStatus::ProcessEvent(AliTRDRawStream *rawStream, Bool_t nocheck)
 251 {
 252   //
 253   // Event Processing loop - AliTRDRawStream
 254   //
 255
 256
 257   Bool_t withInput = kFALSE;
 258
 259   if(!nocheck) {
 260     while (rawStream->Next()) {
 261       Int_t rawversion = rawStream->GetRawVersion();                     //  current raw version
 262       if(rawversion != 2) return kFALSE;
 263       Int_t idetector  = rawStream->GetDet();                            //  current detector
 264       Int_t iRow       = rawStream->GetRow();                            //  current row
 265       Int_t iRowMax    = rawStream->GetMaxRow();                         //  current rowmax
 266       Int_t iCol       = rawStream->GetCol();                            //  current col
 267       Int_t iTimeBin   = rawStream->GetTimeBin();                        //  current time bin
 268       Int_t *signal    = rawStream->GetSignals();                        //  current ADC signal
 269       Int_t nbtimebin  = rawStream->GetNumberOfTimeBins();               //  number of time bins read from data
 270
 271       if((fDetector != -1) && (nbtimebin != fNumberOfTimeBins)) return kFALSE;
 272       fNumberOfTimeBins = nbtimebin;
 273
 274       Int_t fin        = TMath::Min(nbtimebin,(iTimeBin+3));
 275       Int_t n          = 0;
 276
 277       for(Int_t k = iTimeBin; k < fin; k++){
 278         if(signal[n]>0) UpdateHisto(idetector,iRow,iCol,signal[n],iRowMax);
 279         n++;
 280       }
 281
 282       withInput = kTRUE;
 283     }
 284   }
 285   else {
 286     while (rawStream->Next()) {
 287       Int_t idetector  = rawStream->GetDet();                            //  current detector
 288       Int_t iRow       = rawStream->GetRow();                            //  current row
 289       Int_t iRowMax    = rawStream->GetMaxRow();                         //  current rowmax
 290       Int_t iCol       = rawStream->GetCol();                            //  current col
 291       Int_t iTimeBin   = rawStream->GetTimeBin();                        //  current time bin
 292       Int_t *signal    = rawStream->GetSignals();                        //  current ADC signal
 293       Int_t nbtimebin = rawStream->GetNumberOfTimeBins();               //  number of time bins read from data
 294
 295       Int_t fin        = TMath::Min(nbtimebin,(iTimeBin+3));
 296       Int_t n          = 0;
 297
 298       for(Int_t k = iTimeBin; k < fin; k++){
 299         if(signal[n]>0) UpdateHisto(idetector,iRow,iCol,signal[n],iRowMax);
 300         n++;
 301       }
 302
 303       withInput = kTRUE;
 304     }
 305   }
 306
 307   return withInput;
 308 }
 309 //_____________________________________________________________________
 310 Bool_t AliTRDCalibPadStatus::ProcessEvent(AliRawReader *rawReader, Bool_t nocheck)
 311 {
 312   //
 313   //  Event processing loop - AliRawReader
 314   //
 315
 316
 317   AliTRDRawStream rawStream(rawReader);
 318
 319   rawReader->Select("TRD");
 320
 321   return ProcessEvent(&rawStream, nocheck);
 322 }
 323 //_________________________________________________________________________
 324 Bool_t AliTRDCalibPadStatus::ProcessEvent(
 325 #ifdef ALI_DATE
 326                                           eventHeaderStruct *event,
 327                                           Bool_t nocheck
 328 #else
 329                                           eventHeaderStruct* /*event*/,
 330                                           Bool_t /*nocheck*/
 331
 332 #endif
 333                                           )
 334 {
 335   //
 336   //  process date event
 337   //
 338 #ifdef ALI_DATE
 339     AliRawReader *rawReader = new AliRawReaderDate((void*)event);
 340     Bool_t result=ProcessEvent(rawReader, nocheck);
 341     delete rawReader;
 342     return result;
 343 #else
 344     Fatal("AliTRDCalibPadStatus", "this class was compiled without DATE");
 345     return 0;
 346 #endif
 347
 348 }
 349 //_____________________________________________________________________
 350 Bool_t AliTRDCalibPadStatus::TestEvent(Int_t nevent) /*FOLD00*/
 351 {
 352   //
 353   //  Test event loop
 354   // fill one oroc and one iroc with random gaus
 355   //
 356
 357   gRandom->SetSeed(0);
 358
 359     for (Int_t ism=0; ism<18; ism++){
 360         for (Int_t ich=0; ich < 5; ich++){
 361             for (Int_t ipl=0; ipl < 6; ipl++){
 362               for(Int_t irow = 0; irow < fGeo->GetRowMax(ipl,ich,ism); irow++){
 363                 for(Int_t icol = 0; icol < fGeo->GetColMax(ipl); icol++){
 364                   for (Int_t iTimeBin=0; iTimeBin<(30*nevent); iTimeBin++){
 365                     Int_t signal=(Int_t)(ich+8+gRandom->Gaus(0,1.2));
 366                     if ( signal>0 )Update((ipl+ich*6+ism*6*5),irow,icol,signal,fGeo->GetRowMax(ipl,ich,ism));
 367                   }
 368                 }
 369               }
 370             }
 371         }
 372     }
 373     return kTRUE;
 374 }
 375 //_____________________________________________________________________
 376 Bool_t AliTRDCalibPadStatus::TestEventHisto(Int_t nevent) /*FOLD00*/
 377 {
 378   //
 379   //  Test event loop
 380   // fill one oroc and one iroc with random gaus
 381   //
 382
 383   gRandom->SetSeed(0);
 384
 385     for (Int_t ism=0; ism<18; ism++){
 386         for (Int_t ich=0; ich < 5; ich++){
 387             for (Int_t ipl=0; ipl < 6; ipl++){
 388               for(Int_t irow = 0; irow < fGeo->GetRowMax(ipl,ich,ism); irow++){
 389                 for(Int_t icol = 0; icol < fGeo->GetColMax(ipl); icol++){
 390                   for (Int_t iTimeBin=0; iTimeBin<(30*nevent); iTimeBin++){
 391                     Int_t signal=(Int_t)(ich+8+gRandom->Gaus(0,1.2));
 392                     if ( signal>0 )UpdateHisto((ipl+ich*6+ism*6*5),irow,icol,signal,fGeo->GetRowMax(ipl,ich,ism));
 393                   }
 394                 }
 395               }
 396             }
 397         }
 398     }
 399     return kTRUE;
 400 }
 401 //_____________________________________________________________________
 402 TH2F* AliTRDCalibPadStatus::GetHisto(Int_t det, TObjArray *arr, /*FOLD00*/
 403                                   Int_t nbinsY, Float_t ymin, Float_t ymax,
 404                                   Char_t *type, Bool_t force)
 405 {
 406     //
 407     // return pointer to histogram
 408     // if force is true create a new histogram if it doesn't exist allready
 409     //
 410     if ( !force || arr->UncheckedAt(det) )
 411         return (TH2F*)arr->UncheckedAt(det);
 412
 413     // if we are forced and histogram doesn't yes exist create it
 414     Char_t name[255], title[255];
 415
 416     sprintf(name,"hCalib%s%.2d",type,det);
 417     sprintf(title,"%s calibration histogram detector %.2d;ADC channel;Channel (pad)",type,det);
 418
 419     Int_t nbchannels = fGeo->GetRowMax(GetPlane(det),GetChamber(det),GetSector(det))
 420                      * fGeo->GetColMax(GetPlane(det));
 421
 422     // new histogram with calib information. One value for each pad!
 423     TH2F* hist = new TH2F(name,title,
 424                           nbinsY, ymin, ymax,
 425                           nbchannels,0,nbchannels
 426                           );
 427     hist->SetDirectory(0);
 428     arr->AddAt(hist,det);
 429     return hist;
 430 }
 431 //_____________________________________________________________________
 432 TH2F* AliTRDCalibPadStatus::GetHisto(Int_t det, Bool_t force) /*FOLD00*/
 433 {
 434     //
 435     // return pointer to histogram
 436     // if force is true create a new histogram if it doesn't exist allready
 437     //
 438     TObjArray *arr = &fHistoArray;
 439     return GetHisto(det, arr, fAdcMax-fAdcMin, fAdcMin, fAdcMax, "Pedestal", force);
 440 }
 441 //_____________________________________________________________________
 442 AliTRDarrayF* AliTRDCalibPadStatus::GetCal(Int_t det, TObjArray* arr, Bool_t force) /*FOLD00*/
 443 {
 444     //
 445     // return pointer to ROC Calibration
 446     // if force is true create a new AliTRDarrayF if it doesn't exist allready
 447     //
 448     if ( !force || arr->UncheckedAt(det) )
 449         return (AliTRDarrayF*)arr->UncheckedAt(det);
 450
 451     // if we are forced and histogram doesn't yes exist create it
 452     AliTRDarrayF *croc = new AliTRDarrayF();
 453     Int_t nbpad = fGeo->GetRowMax(GetPlane(det),GetChamber(det),GetSector(det))
 454                 * fGeo->GetColMax(GetPlane(det));
 455
 456     // new AliTRDCalROC. One value for each pad!
 457     croc->Expand(nbpad);
 458     for(Int_t k = 0; k < nbpad; k++){
 459       croc->AddAt(0.0,k);
 460     }
 461     arr->AddAt(croc,det);
 462     return croc;
 463 }
 464 //_____________________________________________________________________
 465 AliTRDCalROC* AliTRDCalibPadStatus::GetCalRoc(Int_t det, TObjArray* arr, Bool_t force) /*FOLD00*/
 466 {
 467     //
 468     // return pointer to ROC Calibration
 469     // if force is true create a new AliTRDCalROC if it doesn't exist allready
 470     //
 471     if ( !force || arr->UncheckedAt(det) )
 472         return (AliTRDCalROC*)arr->UncheckedAt(det);
 473
 474     // if we are forced and histogram doesn't yes exist create it
 475
 476     // new AliTRDCalROC. One value for each pad!
 477     AliTRDCalROC *croc = new AliTRDCalROC(GetPlane(det),GetChamber(det));
 478     arr->AddAt(croc,det);
 479     return croc;
 480 }
 481 //_____________________________________________________________________
 482 AliTRDarrayF* AliTRDCalibPadStatus::GetCalEntries(Int_t det, Bool_t force) /*FOLD00*/
 483 {
 484     //
 485     // return pointer to Carge ROC Calibration
 486     // if force is true create a new histogram if it doesn't exist allready
 487     //
 488     TObjArray *arr = &fCalArrayEntries;
 489     return GetCal(det, arr, force);
 490 }
 491 //_____________________________________________________________________
 492 AliTRDarrayF* AliTRDCalibPadStatus::GetCalMean(Int_t det, Bool_t force) /*FOLD00*/
 493 {
 494     //
 495     // return pointer to Carge ROC Calibration
 496     // if force is true create a new histogram if it doesn't exist allready
 497     //
 498     TObjArray *arr = &fCalArrayMean;
 499     return GetCal(det, arr, force);
 500 }
 501 //_____________________________________________________________________
 502 AliTRDarrayF* AliTRDCalibPadStatus::GetCalSquares(Int_t det, Bool_t force) /*FOLD00*/
 503 {
 504     //
 505     // return pointer to Carge ROC Calibration
 506     // if force is true create a new histogram if it doesn't exist allready
 507     //
 508     TObjArray *arr = &fCalArraySquares;
 509     return GetCal(det, arr, force);
 510 }
 511 //_____________________________________________________________________
 512 AliTRDCalROC* AliTRDCalibPadStatus::GetCalRocMean(Int_t det, Bool_t force) /*FOLD00*/
 513 {
 514     //
 515     // return pointer to Carge ROC Calibration
 516     // if force is true create a new histogram if it doesn't exist allready
 517     //
 518     TObjArray *arr = &fCalRocArrayMean;
 519     return GetCalRoc(det, arr, force);
 520 }
 521 //_____________________________________________________________________
 522 AliTRDCalROC* AliTRDCalibPadStatus::GetCalRocRMS(Int_t det, Bool_t force) /*FOLD00*/
 523 {
 524     //
 525     // return pointer to Carge ROC Calibration
 526     // if force is true create a new histogram if it doesn't exist allready
 527     //
 528     TObjArray *arr = &fCalRocArrayRMS;
 529     return GetCalRoc(det, arr, force);
 530 }
 531 //_________________________________________________________________________
 532 void AliTRDCalibPadStatus::Analyse() /*FOLD00*/
 533 {
 534   //
 535   // Calcul the rms properly
 536   //
 537
 538   for(Int_t idet = 0; idet < 540; idet++){
 539
 540     // Take the stuff
 541     fCalEntries                 = ((AliTRDarrayF *)GetCalEntries(idet));
 542     fCalMean                    = ((AliTRDarrayF *)GetCalMean(idet));
 543     fCalSquares                 = ((AliTRDarrayF *)GetCalSquares(idet));
 544
 545     if(!fCalEntries) continue;
 546
 547     AliTRDCalROC *calRocMean    = ((AliTRDCalROC *)GetCalRocMean(idet,kTRUE));
 548     AliTRDCalROC *calRocRMS     = ((AliTRDCalROC *)GetCalRocRMS(idet,kTRUE));
 549
 550     // range channels
 551     Int_t channels = calRocMean->GetNchannels();
 552
 553     for(Int_t ichannels = 0 ; ichannels < channels; ichannels++){
 554
 555       Float_t entries  = fCalEntries->At(ichannels);
 556       Float_t mean     = fCalMean->At(ichannels);
 557       Float_t squares  = fCalSquares->At(ichannels);
 558
 559       Float_t rms = 0.0;
 560       if(entries > 0){
 561         Double_t rm = TMath::Abs(squares-(mean*mean));
 562         rms = TMath::Sqrt(rm);
 563         calRocRMS->SetValue(ichannels,rms/10.0);
 564         calRocMean->SetValue(ichannels,mean/10.0);
 565       }
 566     }
 567   }
 568 }
 569 //_____________________________________________________________________
 570 void AliTRDCalibPadStatus::AnalyseHisto() /*FOLD00*/
 571 {
 572     //
 573     //  Calculate calibration constants
 574     //
 575
 576     Int_t nbinsAdc = fAdcMax-fAdcMin;
 577
 578     TVectorD param(3);
 579     TMatrixD dummy(3,3);
 580
 581     Float_t *array_hP=0;
 582
 583
 584     for (Int_t idet=0; idet<540; idet++){
 585         TH2F *hP = GetHisto(idet);
 586         if ( !hP ) {
 587           continue;
 588         }
 589
 590         AliTRDCalROC *rocMean     = GetCalRocMean(idet,kTRUE);
 591         AliTRDCalROC *rocRMS      = GetCalRocRMS(idet,kTRUE);
 592
 593         array_hP = hP->GetArray();
 594         Int_t nChannels = rocMean->GetNchannels();
 595
 596         for (Int_t iChannel=0; iChannel<nChannels; iChannel++){
 597             Int_t offset = (nbinsAdc+2)*(iChannel+1)+1;
 598             Double_t ret = AliMathBase::FitGaus(array_hP+offset,nbinsAdc,fAdcMin,fAdcMax,&param,&dummy);
 599             // if the fitting failed set noise and pedestal to 0
 600             if ((ret == -4) || (ret == -1)) {
 601                 param[1]=0.0;
 602                 param[2]=0.0;
 603             }
 604             if((param[1]/10.0) > 65534.0) param[1] = 0.0;
 605             if((param[2]/10.0) > 65534.0) param[2] = 0.0;
 606             rocMean->SetValue(iChannel,param[1]/10.0);
 607             rocRMS->SetValue(iChannel,param[2]/10.0);
 608         }
 609     }
 610
 611 }
 612 //_______________________________________________________________________________________
 613 AliTRDCalPadStatus* AliTRDCalibPadStatus::CreateCalPadStatus()
 614 {
 615   //
 616   //
 617   //
 618
 619   AliTRDCalPadStatus* obj = new AliTRDCalPadStatus("padstatus", "padstatus");
 620
 621   for (Int_t idet=0; idet<540; ++idet)
 622     {
 623       AliTRDCalSingleChamberStatus *calROC = obj->GetCalROC(idet);
 624
 625       //Take the stuff
 626       fCalEntries = ((AliTRDarrayF *)GetCalEntries(idet));
 627       AliTRDCalROC *calRocMean    = ((AliTRDCalROC *)GetCalRocMean(idet));
 628       AliTRDCalROC *calRocRMS     = ((AliTRDCalROC *)GetCalRocRMS(idet));
 629
 630       if ( !calRocMean ) {
 631         for(Int_t k = 0; k < calROC->GetNchannels(); k++){
 632           calROC->SetStatus(k,AliTRDCalPadStatus::kMasked);
 633         }
 634         //printf("no fCalRocMean for %d\n",idet);
 635         continue;
 636       }
 637
 638       //Range
 639       Int_t channels = calROC->GetNchannels();
 640
 641
 642       //Mean
 643       Float_t meanentries = 0.0;
 644       if(!fCalEntries){
 645         if(GetHisto(idet)){
 646           meanentries = GetHisto(idet)->GetEntries()/(channels);
 647         }
 648       }
 649       else meanentries       = TMath::Mean(channels,((TArrayF *)fCalEntries)->GetArray());
 650       //Double_t meanmean      = calRocMean->GetMean()*10.0;
 651       //Double_t meansquares   = calRocRMS->GetMean()*10.0;
 652
 653
 654       for(Int_t ich = 0; ich < channels; ich++){
 655
 656         Float_t entries = 0.0;
 657         if(!fCalEntries){
 658           if(GetHisto(idet)){
 659             for(Int_t bin = 0; bin < (fAdcMax-fAdcMin); bin++){
 660               entries += GetHisto(idet)->GetArray()[(ich+1)*(fAdcMax-fAdcMin+2)+(bin+1)];
 661             }
 662           }
 663         }
 664         else entries     = fCalEntries->At(ich);
 665         //Float_t mean     = calRocMean->GetValue(ich)*10.0;
 666         Float_t rms      = calRocRMS->GetValue(ich)*10.0;
 667         if(ich > 1720) printf("rms %f\n",rms);
 668
 669         //if((entries < 0.3*meanentries) || (TMath::Abs(rms-meansquares) > ((Float_t)(fAdcMax-fAdcMin)/2.0)) || (rms == 0.0)) calROC->SetStatus(ich, AliTRDCalPadStatus::kMasked);
 670         if(rms <= 0.01) calROC->SetStatus(ich, AliTRDCalPadStatus::kMasked);
 671         }
 672     }
 673
 674   return obj;
 675
 676 }
 677 //_____________________________________________________________________
 678 void AliTRDCalibPadStatus::DumpToFile(const Char_t *filename, const Char_t *dir, Bool_t append) /*FOLD00*/
 679 {
 680     //
 681     //  Write class to file
 682     //
 683
 684     TString sDir(dir);
 685     TString option;
 686
 687     if ( append )
 688         option = "update";
 689     else
 690         option = "recreate";
 691
 692     TDirectory *backup = gDirectory;
 693     TFile f(filename,option.Data());
 694     f.cd();
 695     if ( !sDir.IsNull() ){
 696         f.mkdir(sDir.Data());
 697         f.cd(sDir);
 698     }
 699     this->Write();
 700     f.Close();
 701
 702     if ( backup ) backup->cd();
 703 }//_____________________________________________________________________________
 704 Int_t AliTRDCalibPadStatus::GetPlane(Int_t d) const
 705 {
 706   //
 707   // Reconstruct the plane number from the detector number
 708   //
 709
 710   return ((Int_t) (d % 6));
 711
 712 }
 713
 714 //_____________________________________________________________________________
 715 Int_t AliTRDCalibPadStatus::GetChamber(Int_t d) const
 716 {
 717   //
 718   // Reconstruct the chamber number from the detector number
 719   //
 720
 721   return ((Int_t) (d % 30) / 6);
 722
 723 }
 724 //_____________________________________________________________________________
 725 Int_t AliTRDCalibPadStatus::GetSector(Int_t d) const
 726 {
 727   //
 728   // Reconstruct the sector number from the detector number
 729   //
 730
 731   return ((Int_t) (d / 30));
 732
 733 }
 734
 735