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 //                                                                        //
  20 // Example: fill pedestal with Gaussian noise                             //
  21 //                                                                        //
  22 //  AliTRDCalibPadStatus ped;                                             //
  23 //  ped.TestEvent(numberofevent);                                         //
  24 //                                                                        //
  25 //  // Method without histo                                               //
  26 //  ped.Analyse();                                                        //
  27 //                                                                        //
  28 //  // Create the histo of the AliTRDCalROC                               //
  29 //  TH2F * histo2dm = ped.GetCalRocMean(0,kFALSE)->MakeHisto2D();         //
  30 //  histo2dm->Scale(10.0);                                                //
  31 //  TH1F * histo1dm = ped.GetCalRocMean(0,kFALSE)->MakeHisto1D();         //
  32 //  histo1dm->Scale(10.0);                                                //
  33 //  TH2F * histo2ds = ped.GetCalRocSquares(0,kFALSE)->MakeHisto2D();      //
  34 //  histo2ds->Scale(10.0);                                                //
  35 //  TH1F * histo1ds = ped.GetCalRocSquares(0,kFALSE)->MakeHisto1D();      //
  36 //  histo1ds->Scale(10.0)                                                 //
  37 //                                                                        //
  38 //  // Draw output                                                        //
  39 //  TCanvas* c1 = new TCanvas;                                            //
  40 //  c1->Divide(2,2);                                                      //
  41 //  c1->cd(1);                                                            //
  42 //  histo2dm->Draw("colz");                                               //
  43 //  c1->cd(2);                                                            //
  44 //  histo1dm->Draw();                                                     //
  45 //  c1->cd(3);                                                            //
  46 //  histo2ds->Draw("colz");                                               //
  47 //  c1->cd(4);                                                            //
  48 //  histo1ds->Draw();                                                     //
  49 //                                                                        //
  50 //  // Method with histo                                                  //
  51 //  ped.AnalyseHisto();                                                   //
  52 //                                                                        //
  53 //  // Take the histo                                                     //
  54 //  TH1F *histo = ped.GetHisto(31);                                       //
  55 //  histo->SetEntries(1);                                                 //
  56 //  histo->Draw();                                                        //
  57 //
  58 // Authors:
  59 //   R. Bailhache (R.Bailhache@gsi.de, rbailhache@ikf.uni-frankfurt.de)
  60 //   J. Book (jbook@ikf.uni-frankfurt.de)
  61 //                                                                                                    //
  62 ////////////////////////////////////////////////////////////////////////////
  63
  64
  65 //Root includes
  66 #include <TObjArray.h>
  67 #include <TH2F.h>
  68 #include <TString.h>
  69 #include <TMath.h>
  70 #include <TRandom.h>
  71
  72 //#include <TRandom.h>
  73 #include <TDirectory.h>
  74 #include <TFile.h>
  75
  76 //AliRoot includes
  77 #include <AliMathBase.h>
  78 #include "AliRawReader.h"
  79 #include "AliRawReaderRoot.h"
  80 #include "AliRawReaderDate.h"
  81
  82 //header file
  83 #include "AliLog.h"
  84 #include "AliTRDCalibPadStatus.h"
  85 #include "AliTRDrawStreamBase.h"
  86 #include "AliTRDgeometry.h"
  87 #include "AliTRDCommonParam.h"
  88 #include "./Cal/AliTRDCalROC.h"
  89 #include "./Cal/AliTRDCalPadStatus.h"
  90 #include "./Cal/AliTRDCalDet.h"
  91 #include "./Cal/AliTRDCalPad.h"
  92 #include "./Cal/AliTRDCalSingleChamberStatus.h"
  93
  94 #include "AliTRDrawFastStream.h"
  95 #include "AliTRDdigitsManager.h"
  96 #include "AliTRDdigitsParam.h"
  97 #include "AliTRDSignalIndex.h"
  98 #include "AliTRDarraySignal.h"
  99 #include "AliTRDarrayADC.h"
 100 #include "AliTRDfeeParam.h"
 101
 102 #ifdef ALI_DATE
 103 #include "event.h"
 104 #endif
 105
 106 ClassImp(AliTRDCalibPadStatus) /*FOLD00*/
 107
 108 //_____________________________________________________________________
 109 AliTRDCalibPadStatus::AliTRDCalibPadStatus() : /*FOLD00*/
 110   TObject(),
 111   fGeo(0),
 112   fAdcMin(0),
 113   fAdcMax(21),
 114   fDetector(-1),
 115   fNumberOfTimeBins(0),
 116   fCalRocArrayMean(540),
 117   fCalRocArrayRMS(540),
 118   fCalRocArrayMeand(540),
 119   fCalRocArrayRMSd(540),
 120   fHistoArray(540)
 121 {
 122     //
 123     // default constructor
 124     //
 125
 126   fGeo = new AliTRDgeometry();
 127
 128 }
 129
 130 //_____________________________________________________________________
 131 AliTRDCalibPadStatus::AliTRDCalibPadStatus(const AliTRDCalibPadStatus &ped) : /*FOLD00*/
 132   TObject(ped),
 133   fGeo(0),
 134   fAdcMin(ped.GetAdcMin()),
 135   fAdcMax(ped.GetAdcMax()),
 136   fDetector(ped.fDetector),
 137   fNumberOfTimeBins(ped.fNumberOfTimeBins),
 138   fCalRocArrayMean(540),
 139   fCalRocArrayRMS(540),
 140   fCalRocArrayMeand(540),
 141   fCalRocArrayRMSd(540),
 142   fHistoArray(540)
 143 {
 144     //
 145     // copy constructor
 146     //
 147     for (Int_t idet = 0; idet < 540; idet++){
 148         const AliTRDCalROC *calRocMean  = (AliTRDCalROC*)ped.fCalRocArrayMean.UncheckedAt(idet);
 149         const AliTRDCalROC *calRocRMS   = (AliTRDCalROC*)ped.fCalRocArrayRMS.UncheckedAt(idet);
 150         const AliTRDCalROC *calRocMeand = (AliTRDCalROC*)ped.fCalRocArrayMeand.UncheckedAt(idet);
 151         const AliTRDCalROC *calRocRMSd  = (AliTRDCalROC*)ped.fCalRocArrayRMSd.UncheckedAt(idet);
 152         const TH2F         *hped        = (TH2F*)ped.fHistoArray.UncheckedAt(idet);
 153
 154         if ( calRocMean != 0x0 ) fCalRocArrayMean.AddAt(new AliTRDCalROC(*calRocMean), idet);
 155         if ( calRocRMS != 0x0 )  fCalRocArrayRMS.AddAt(new AliTRDCalROC(*calRocRMS), idet);
 156
 157         if ( calRocMeand != 0x0 ) fCalRocArrayMeand.AddAt(new AliTRDCalROC(*calRocMeand), idet);
 158         if ( calRocRMSd != 0x0 )  fCalRocArrayRMSd.AddAt(new AliTRDCalROC(*calRocRMSd), idet);
 159
 160         if ( hped != 0x0 ){
 161           TH2F *hNew = new TH2F(*hped);
 162           hNew->SetDirectory(0);
 163           fHistoArray.AddAt(hNew,idet);
 164         }
 165
 166     }
 167     if (fGeo) {
 168       delete fGeo;
 169     }
 170     fGeo = new AliTRDgeometry();
 171 }
 172
 173 //_____________________________________________________________________
 174 AliTRDCalibPadStatus& AliTRDCalibPadStatus::operator = (const  AliTRDCalibPadStatus &source)
 175 {
 176   //
 177   // assignment operator
 178   //
 179   if (&source == this) return *this;
 180   new (this) AliTRDCalibPadStatus(source);
 181
 182   return *this;
 183 }
 184 //_____________________________________________________________________
 185 AliTRDCalibPadStatus::~AliTRDCalibPadStatus() /*FOLD00*/
 186 {
 187   //
 188   // destructor
 189   //
 190   fCalRocArrayMean.Delete();
 191   fCalRocArrayRMS.Delete();
 192   fCalRocArrayMeand.Delete();
 193   fCalRocArrayRMSd.Delete();
 194   fHistoArray.Delete();
 195   if (fGeo) {
 196     delete fGeo;
 197   }
 198 }
 199 //_____________________________________________________________________
 200 void AliTRDCalibPadStatus::Destroy()
 201 {
 202   //
 203   // Destroy
 204   //
 205   fCalRocArrayMean.Delete();
 206   fCalRocArrayRMS.Delete();
 207   fCalRocArrayMeand.Delete();
 208   fCalRocArrayRMSd.Delete();
 209   fHistoArray.Delete();
 210 }
 211 //_____________________________________________________________________
 212 Int_t AliTRDCalibPadStatus::UpdateHisto(const Int_t icdet, /*FOLD00*/
 213                                         const Int_t icRow,
 214                                         const Int_t icCol,
 215                                         const Int_t csignal,
 216                                         const Int_t crowMax,
 217                                         const Int_t ccold,
 218                                         const Int_t icMcm)
 219 {
 220   //
 221   // Signal filling methode
 222   //
 223   Int_t nbchannel = icRow+icCol*crowMax;
 224
 225   // now the case of double read channel
 226   if(ccold > 0){
 227     nbchannel = (((ccold-1)*8+ icMcm)*crowMax+icRow)+144*crowMax;
 228     //printf("nbchannel %d, ccold %d, icMcm %d, crowMax %d, icRow %d\n",nbchannel,ccold,icMcm,crowMax,icRow);
 229   }
 230
 231   // fast filling methode.
 232   // Attention: the entry counter of the histogram is not increased
 233   //            this means that e.g. the colz draw option gives an empty plot
 234
 235   Int_t bin = 0;
 236
 237   if ( !(((Int_t)csignal>=fAdcMax ) || ((Int_t)csignal<fAdcMin)) )
 238     bin = (nbchannel+1)*(fAdcMax-fAdcMin+2)+((Int_t)csignal-fAdcMin+1);
 239
 240   //GetHisto(icdet,kTRUE)->Fill(csignal,nbchannel);
 241
 242   GetHisto(icdet,kTRUE)->GetArray()[bin]++;
 243
 244   return 0;
 245 }
 246 //_____________________________________________________________________
 247 Int_t AliTRDCalibPadStatus::UpdateHisto2(const Int_t icdet, /*FOLD00*/
 248                                          const Int_t icRow,
 249                                          const Int_t icCol,
 250                                          const Int_t csignal,
 251                                          const Int_t crowMax,
 252                                          const Int_t ccold,
 253                                          const Int_t icMcm,
 254                                          const Int_t icRob
 255                                          )
 256 {
 257   //
 258   // Signal filling methode
 259   //
 260   Int_t nbchannel = icRow+icCol*crowMax;
 261   Int_t mCm = icMcm%4;
 262   Int_t rOb = icRob%2;
 263
 264   // now the case of double read channel
 265   if(ccold > 0){
 266     nbchannel = (((ccold-1)*8+ (mCm+rOb*4))*crowMax+icRow)+144*crowMax;
 267     //printf("nbchannel %d, ccold %d, icMcm %d, crowMax %d, icRow %d\n",nbchannel,ccold,icMcm,crowMax,icRow);
 268   }
 269
 270   // fast filling methode.
 271   // Attention: the entry counter of the histogram is not increased
 272   //            this means that e.g. the colz draw option gives an empty plot
 273
 274   Int_t bin = 0;
 275
 276   if ( !(((Int_t)csignal>=fAdcMax ) || ((Int_t)csignal<fAdcMin)) )
 277     bin = (nbchannel+1)*(fAdcMax-fAdcMin+2)+((Int_t)csignal-fAdcMin+1);
 278
 279   //GetHisto(icdet,kTRUE)->Fill(csignal,nbchannel);
 280
 281   GetHisto(icdet,kTRUE)->GetArray()[bin]++;
 282
 283   return 0;
 284 }
 285 //_____________________________________________________________________
 286 Int_t AliTRDCalibPadStatus::ProcessEvent(AliTRDrawStreamBase *rawStream, Bool_t nocheck)
 287 {
 288   //
 289   // Event Processing loop - AliTRDRawStreamCosmic
 290   // 0 time bin problem or zero suppression
 291   // 1 no input
 292   // 2 input
 293   //
 294
 295   //
 296   // Raw version number:
 297   // [3,31] non zero suppressed
 298   // 2,4 and [32,63] zero suppressed
 299   //
 300
 301   Int_t withInput = 1;
 302
 303   rawStream->SetSharedPadReadout(kTRUE);
 304
 305   if(!nocheck) {
 306
 307     // Check the raw version and if all have the same number of timebins.
 308
 309     while (rawStream->Next()) {
 310
 311       Int_t rawversion = rawStream->GetRawVersion();                     //  current raw version
 312       //printf("Raw version is %d\n",rawversion);
 313
 314       // Could eventually change, have to check with time
 315       if((rawversion < 3) || (rawversion > 31)) {
 316         AliInfo(Form("this is not no-zero-suppressed data, the version is %d",rawversion));
 317         return 0;
 318       }
 319       Int_t idetector  = rawStream->GetDet();                            //  current detector
 320       Int_t iRow       = rawStream->GetRow();                            //  current row
 321       Int_t iRowMax    = rawStream->GetMaxRow();                         //  current rowmax
 322       Int_t iCol       = rawStream->GetCol();                            //  current col
 323       Int_t iADC       = 21-rawStream->GetADC();                         //  current ADC
 324
 325       // It goes in the opposite direction
 326       Int_t col        = 0;
 327       if(iADC == 1) col = 1;
 328       else {
 329         col = TMath::Max(0,(Int_t)(iADC-19));
 330         if(col > 0) col++;
 331       }
 332       Int_t mcm        = (Int_t)(iCol/18);                               //  current group of 18 col pads
 333       if(col > 1) mcm -= 1;
 334       if(col ==1) mcm += 1;
 335
 336       // printf to check
 337       //Bool_t shared = rawStream->IsCurrentPadShared();
 338       //printf("ADC %d, iCol %d, col %d, mcm %d, shared %d\n",iADC,iCol,col,mcm,(Int_t)shared);
 339
 340       // Take the signal
 341       Int_t *signal    = rawStream->GetSignals();                        //  current ADC signal
 342       Int_t nbtimebin  = rawStream->GetNumberOfTimeBins();               //  number of time bins read from data
 343
 344       if((fDetector != -1) && (nbtimebin != fNumberOfTimeBins)) {
 345         AliInfo(Form("the number of time bins is %d, is different from the previous one %d",nbtimebin,fNumberOfTimeBins));
 346         return 0;
 347       }
 348       fNumberOfTimeBins = nbtimebin;
 349       fDetector         = idetector;
 350
 351       for(Int_t k = 0; k < fNumberOfTimeBins; k++){
 352         if(signal[k]>0 && iCol != -1) UpdateHisto(idetector,iRow,iCol,signal[k],iRowMax,col,mcm);
 353       }
 354
 355       withInput = 2;
 356     }
 357   }
 358   else {
 359
 360     while (rawStream->Next()) {
 361
 362       Int_t idetector  = rawStream->GetDet();                            //  current detector
 363       Int_t iRow       = rawStream->GetRow();                            //  current row
 364       Int_t iRowMax    = rawStream->GetMaxRow();                         //  current rowmax
 365       Int_t iCol       = rawStream->GetCol();                            //  current col
 366       Int_t iADC       = 21-rawStream->GetADC();                            //  current ADC
 367
 368       // It goes in the opposite direction
 369       Int_t col        = 0;
 370       if(iADC == 1) col = 1;
 371       else {
 372         col = TMath::Max(0,(Int_t)(iADC-19));
 373         if(col > 0) col++;
 374       }
 375       Int_t mcm        = (Int_t)(iCol/18);                               //  current group of 18 col pads
 376       if(col > 1) mcm -= 1;
 377       if(col ==1) mcm += 1;
 378
 379       // Take the signal
 380       Int_t *signal    = rawStream->GetSignals();                        //  current ADC signal
 381       Int_t nbtimebin = rawStream->GetNumberOfTimeBins();               //  number of time bins read from data
 382
 383
 384       //printf("det %d, row %d, signal[0] %d, signal[1] %d, signal [2] %d\n", idetector, iRow, signal[0], signal[1], signal[2]);
 385
 386       for(Int_t k = 0; k < nbtimebin; k++){
 387         if(signal[k]>0 && iCol != -1) {
 388           UpdateHisto(idetector,iRow,iCol,signal[k],iRowMax,col,mcm);
 389           //printf("Update with det %d, row %d, col %d, signal %d, rowmax %d, col %d, mcm %d\n",idetector,iRow,iCol,signal[n],iRowMax,col,mcm);
 390         }
 391       }
 392
 393       withInput = 2;
 394     }
 395   }
 396
 397   return withInput;
 398 }
 399 //_____________________________________________________________________
 400 Int_t AliTRDCalibPadStatus::ProcessEvent(AliRawReader *rawReader, Bool_t nocheck)
 401 {
 402   //
 403   //  Event processing loop - AliRawReader
 404   //
 405
 406   Int_t result;
 407
 408   rawReader->Select("TRD");
 409
 410   AliTRDrawStreamBase *pstream = AliTRDrawStreamBase::GetRawStream(rawReader);
 411
 412   result = ProcessEvent(pstream, nocheck);
 413
 414   delete pstream;
 415
 416   return result;
 417 }
 418
 419 //_________________________________________________________________________
 420 Int_t AliTRDCalibPadStatus::ProcessEvent(
 421 #ifdef ALI_DATE
 422                                           const eventHeaderStruct *event,
 423                                           Bool_t nocheck
 424 #else
 425                                           const eventHeaderStruct* /*event*/,
 426                                           Bool_t /*nocheck*/
 427
 428 #endif
 429                                           )
 430 {
 431   //
 432   //  process date event
 433   //
 434 #ifdef ALI_DATE
 435     AliRawReader *rawReader = new AliRawReaderDate((void*)event);
 436     Bool_t result=ProcessEvent(rawReader, nocheck);
 437     delete rawReader;
 438     return result;
 439 #else
 440     Fatal("AliTRDCalibPadStatus", "this class was compiled without DATE");
 441     return 0;
 442 #endif
 443
 444 }
 445
 446 //_____________________________________________________________________
 447 Int_t AliTRDCalibPadStatus::ProcessEvent2(AliRawReader *rawReader)
 448 {
 449   //
 450   // Event Processing loop - AliTRDRawStreamCosmic
 451   // 0 time bin problem or zero suppression
 452   // 1 no input
 453   // 2 input
 454   // Raw version number:
 455   // [3,31] non zero suppressed
 456   // 2,4 and [32,63] zero suppressed
 457   //
 458
 459   Int_t withInput = 1;
 460
 461   AliTRDrawFastStream *rawStream = new AliTRDrawFastStream(rawReader);
 462   rawStream->SetSharedPadReadout(kTRUE);
 463
 464   AliTRDdigitsManager *digitsManager = new AliTRDdigitsManager(kTRUE);
 465   digitsManager->CreateArrays();
 466
 467   AliTRDfeeParam *feeParam = AliTRDfeeParam::Instance();
 468
 469   Int_t det    = 0;
 470   while ((det = rawStream->NextChamber(digitsManager, NULL, NULL)) >= 0) { //idetector
 471     if (digitsManager->GetIndexes(det)->HasEntry()) {//QA
 472       //        printf("there is ADC data on this chamber!\n");
 473
 474        AliTRDarrayADC *digits = (AliTRDarrayADC *) digitsManager->GetDigits(det); //mod
 475       if (digits->HasData()) { //array
 476
 477         AliTRDSignalIndex   *indexes = digitsManager->GetIndexes(det);
 478         if (indexes->IsAllocated() == kFALSE) {
 479           AliError("Indexes do not exist!");
 480           break;
 481         }
 482         Int_t iRow  = 0;
 483         Int_t iCol  = 0;
 484         indexes->ResetCounters();
 485
 486         while (indexes->NextRCIndex(iRow, iCol)) { //column,row
 487
 488
 489           AliTRDdigitsParam *digitParam = (AliTRDdigitsParam *)digitsManager->GetDigitsParam();
 490
 491           Int_t mcm          = 0;     // MCM from AliTRDfeeParam
 492           Int_t rob          = 0;     // ROB from AliTRDfeeParam
 493           Int_t extCol       = 0;     // extended column from  AliTRDfeeParam
 494           mcm = feeParam->GetMCMfromPad(iRow,iCol);
 495           rob = feeParam->GetROBfromPad(iRow,iCol);
 496
 497           Int_t idetector  = det;                            //  current detector
 498           Int_t iRowMax    = rawStream->GetMaxRow();         //  current rowmax
 499
 500           Int_t adc        = 20 - (iCol%18) -1;                 //  current adc
 501           Int_t col        = 0;                              //  col!=0 ->Shared Pad
 502           extCol = feeParam->GetExtendedPadColFromADC(rob,mcm,adc);
 503           //printf("  iCol %d  iRow %d  iRowMax %d  rob %d  mcm %d  adc %d  extCol %d\n",iCol,iRow,iRowMax,rob,mcm,adc,extCol);
 504
 505           // Signal for regular pads
 506           Int_t nbtimebin  = digitParam->GetNTimeBins();     //  number of time bins read from data
 507           for(Int_t k = 0; k < nbtimebin; k++){
 508             Short_t signal = 0;
 509             signal = digits->GetData(iRow,iCol,k);
 510
 511             if(signal>0) {
 512               UpdateHisto2(idetector,iRow,iCol,signal,iRowMax,col,mcm,rob);
 513             }
 514           }
 515
 516
 517
 518           if((adc==3-1 || adc==20-1 || adc==19-1) && (iCol > 1 && iCol <142) /* && fSharedPadsOn*/ ) { //SHARED PADS
 519
 520               switch(adc) {
 521               case 2:
 522                 adc = 20;                                       //shared Pad adc
 523                 mcm = feeParam->GetMCMfromSharedPad(iRow,iCol); //shared Pad mcm
 524                 col =  1;
 525                 break;
 526               case 19:
 527                 adc = 1;                                        //shared Pad adc
 528                 mcm = feeParam->GetMCMfromSharedPad(iRow,iCol); //shared Pad mcm
 529                 col =  2;
 530                 break;
 531               case 18:
 532                 adc =  0;                                       //shared Pad adc
 533                 mcm = feeParam->GetMCMfromSharedPad(iRow,iCol); //shared Pad mcm
 534                 col =  3;
 535                 break;
 536               }
 537               rob = feeParam->GetROBfromSharedPad(iRow,iCol);     //shared Pad rob
 538
 539
 540             extCol = feeParam->GetExtendedPadColFromADC(rob,mcm,adc);     //extended pad col via the shared pad rob,mcm and adc
 541
 542             //printf("SHARED PAD ---  iCol %d  iRow %d  rob %d  mcm %d  adc %d  extCol %d  col %d\n",iCol,iRow,rob,mcm,adc,extCol,col);
 543             for(Int_t k = 0; k < nbtimebin; k++){
 544               Short_t signal = 0;
 545               signal = digits->GetDataByAdcCol(iRow,extCol,k);
 546
 547               if(signal>0) {
 548                 UpdateHisto2(idetector,iRow,iCol,signal,iRowMax,col,mcm,rob);
 549               }
 550             }
 551           } //shared pads end
 552
 553
 554           withInput = 2;
 555         }//column,row
 556
 557       }//array
 558     }//QA
 559     digitsManager->ClearArrays(det);
 560   }//idetector
 561   delete digitsManager;
 562   delete rawStream;
 563   return withInput;
 564 }
 565
 566 //_____________________________________________________________________
 567 Bool_t AliTRDCalibPadStatus::TestEventHisto(Int_t nevent, Int_t sm, Int_t ch) /*FOLD00*/
 568 {
 569   //
 570   //  Test event loop
 571   // fill one oroc and one iroc with random gaus
 572   //
 573
 574   gRandom->SetSeed(0);
 575
 576     for (Int_t ism=sm; ism<sm+1; ism++){
 577         for (Int_t ich=ch; ich < ch+1; ich++){
 578             for (Int_t ipl=0; ipl < 6; ipl++){
 579               for(Int_t irow = 0; irow < fGeo->GetRowMax(ipl,ich,ism); irow++){
 580                 for(Int_t icol = 0; icol < fGeo->GetColMax(ipl); icol++){
 581                   for (Int_t iTimeBin=0; iTimeBin<(30*nevent); iTimeBin++){
 582                     Int_t signal=TMath::Nint(gRandom->Gaus(10,1.5));
 583                     if ( signal>0 )UpdateHisto((ipl+ich*6+ism*6*5),irow,icol,signal,fGeo->GetRowMax(ipl,ich,ism),0,0);
 584                   }
 585                 }
 586               }
 587             }
 588         }
 589     }
 590     return kTRUE;
 591 }
 592
 593 //_____________________________________________________________________
 594 TH2F* AliTRDCalibPadStatus::GetHisto(Int_t det, TObjArray *arr, /*FOLD00*/
 595                                   Int_t nbinsY, Float_t ymin, Float_t ymax,
 596                                   const Char_t *type, Bool_t force)
 597 {
 598     //
 599     // return pointer to histogram
 600     // if force is true create a new histogram if it doesn't exist allready
 601     //
 602     if ( !force || arr->UncheckedAt(det) )
 603         return (TH2F*)arr->UncheckedAt(det);
 604
 605     // if we are forced and histogram doesn't yes exist create it
 606     Char_t name[255], title[255];
 607
 608     sprintf(name,"hCalib%s%.3d",type,det);
 609     sprintf(title,"%s calibration histogram detector %.2d;ADC channel;Channel (pad)",type,det);
 610
 611
 612     Int_t nbchannels = fGeo->GetRowMax(GetLayer(det),GetStack(det),GetSector(det))*fGeo->GetColMax(GetLayer(det));
 613
 614     // we will add 3*8*rowMax channels at the end for the double counted
 615     nbchannels += 3*8*(fGeo->GetRowMax(GetLayer(det),GetStack(det),GetSector(det)));
 616
 617
 618     // new histogram with calib information. One value for each pad!
 619     TH2F* hist = new TH2F(name,title,
 620                           nbinsY, ymin, ymax,
 621                           nbchannels,0,nbchannels
 622                           );
 623     hist->SetDirectory(0);
 624     arr->AddAt(hist,det);
 625     return hist;
 626 }
 627
 628 //_____________________________________________________________________
 629 TH2F* AliTRDCalibPadStatus::GetHisto(Int_t det, Bool_t force) /*FOLD00*/
 630 {
 631     //
 632     // return pointer to histogram
 633     // if force is true create a new histogram if it doesn't exist allready
 634     //
 635     TObjArray *arr = &fHistoArray;
 636     return GetHisto(det, arr, fAdcMax-fAdcMin, fAdcMin-0.5, fAdcMax-0.5, "Pedestal", force);
 637 }
 638
 639 //_____________________________________________________________________
 640 AliTRDCalROC* AliTRDCalibPadStatus::GetCalRoc(Int_t det, TObjArray* arr, Bool_t force) /*FOLD00*/
 641 {
 642     //
 643     // return pointer to ROC Calibration
 644     // if force is true create a new AliTRDCalROC if it doesn't exist allready
 645     //
 646     if ( !force || arr->UncheckedAt(det) )
 647         return (AliTRDCalROC*)arr->UncheckedAt(det);
 648
 649     // if we are forced and histogram doesn't yes exist create it
 650
 651     // new AliTRDCalROC. One value for each pad!
 652     AliTRDCalROC *croc = new AliTRDCalROC(GetLayer(det),GetStack(det));
 653     arr->AddAt(croc,det);
 654     return croc;
 655 }
 656 //_____________________________________________________________________
 657 AliTRDCalROC* AliTRDCalibPadStatus::GetCalRocMean(Int_t det, Bool_t force) /*FOLD00*/
 658 {
 659     //
 660     // return pointer to Carge ROC Calibration
 661     // if force is true create a new histogram if it doesn't exist allready
 662     //
 663     TObjArray *arr = &fCalRocArrayMean;
 664     return GetCalRoc(det, arr, force);
 665 }
 666
 667 //_____________________________________________________________________
 668 AliTRDCalROC* AliTRDCalibPadStatus::GetCalRocRMS(Int_t det, Bool_t force) /*FOLD00*/
 669 {
 670     //
 671     // return pointer to Carge ROC Calibration
 672     // if force is true create a new histogram if it doesn't exist allready
 673     //
 674     TObjArray *arr = &fCalRocArrayRMS;
 675     return GetCalRoc(det, arr, force);
 676 }
 677 //_____________________________________________________________________
 678 AliTRDCalROC* AliTRDCalibPadStatus::GetCalRocMeand(Int_t det, Bool_t force) /*FOLD00*/
 679 {
 680     //
 681     // return pointer to Carge ROC Calibration
 682     // if force is true create a new histogram if it doesn't exist allready
 683     //
 684     TObjArray *arr = &fCalRocArrayMeand;
 685     return GetCalRoc(det, arr, force);
 686 }
 687
 688 //_____________________________________________________________________
 689 AliTRDCalROC* AliTRDCalibPadStatus::GetCalRocRMSd(Int_t det, Bool_t force) /*FOLD00*/
 690 {
 691     //
 692     // return pointer to Carge ROC Calibration
 693     // if force is true create a new histogram if it doesn't exist allready
 694     //
 695     TObjArray *arr = &fCalRocArrayRMSd;
 696     return GetCalRoc(det, arr, force);
 697 }
 698
 699 //_____________________________________________________________________
 700 void AliTRDCalibPadStatus::AnalyseHisto() /*FOLD00*/
 701 {
 702     //
 703     //  Calculate calibration constants
 704     //
 705
 706     Int_t nbinsAdc = fAdcMax-fAdcMin;
 707
 708     TVectorD param(4);
 709     TMatrixD dummy(3,3);
 710
 711     Float_t *arrayHP=0;
 712
 713
 714     for (Int_t idet=0; idet<540; idet++){
 715         TH2F *hP = GetHisto(idet);
 716         if ( !hP ) {
 717           continue;
 718         }
 719
 720         //printf("Entries for %d\n",idet);
 721
 722         AliTRDCalROC *rocMean     = GetCalRocMean(idet,kTRUE);
 723         AliTRDCalROC *rocRMS      = GetCalRocRMS(idet,kTRUE);
 724
 725         arrayHP = hP->GetArray();
 726         Int_t nChannels = rocMean->GetNchannels();
 727
 728         for (Int_t iChannel=0; iChannel<nChannels; iChannel++){
 729             Int_t offset = (nbinsAdc+2)*(iChannel+1)+1;
 730             Double_t ret = AliMathBase::FitGaus(arrayHP+offset,nbinsAdc,fAdcMin-0.5,fAdcMax-0.5,&param,&dummy);
 731             // if the fitting failed set noise and pedestal to 0
 732             if ((ret==-4) || (ret==-1) || (ret==-2)) {
 733                 param[1]=0.0;
 734                 param[2]=0.0;
 735             }
 736             if((param[1]/10.0) > 65534.0) param[1] = 0.0;
 737             if((param[2]/10.0) > 65534.0) param[2] = 0.0;
 738             rocMean->SetValue(iChannel,param[1]/10.0);
 739             rocRMS->SetValue(iChannel,param[2]/10.0);
 740         }
 741
 742         // here we analyse doubled read channels
 743
 744         AliTRDCalROC *rocMeand     = GetCalRocMeand(idet,kTRUE);
 745         AliTRDCalROC *rocRMSd      = GetCalRocRMSd(idet,kTRUE);
 746
 747         Int_t nrows = rocMeand->GetNrows();
 748         Int_t shift = 144*nrows;
 749         Int_t total = shift+3*8*nrows;
 750
 751         for (Int_t iChannel=shift; iChannel<total; iChannel++){
 752             Int_t offset = (nbinsAdc+2)*(iChannel+1)+1;
 753             Double_t ret = AliMathBase::FitGaus(arrayHP+offset,nbinsAdc,fAdcMin-0.5,fAdcMax-0.5,&param,&dummy);
 754             // if the fitting failed set noise and pedestal to 0
 755             if ((ret==-4) || (ret==-1) || (ret==-2)) {
 756                 param[1]=0.0;
 757                 param[2]=0.0;
 758             }
 759             if((param[1]/10.0) > 65534.0) param[1] = 0.0;
 760             if((param[2]/10.0) > 65534.0) param[2] = 0.0;
 761
 762             // here we have to recalculate backward
 763             Int_t nb   = iChannel-shift;
 764             Int_t row  = nb%nrows;
 765             Int_t j    = (Int_t)(nb/nrows);
 766             Int_t imcm = j%8;
 767             Int_t icol = (Int_t)(j/8);
 768
 769             Int_t finalcol = 18*imcm;
 770             if(icol > 0) icol += 17;
 771             else icol = -1;
 772             finalcol += icol;
 773
 774             Int_t channel = row+finalcol*nrows;
 775
 776             //printf("iChannel %d, nrows %d, finalcol %d, row %d, channel %d\n",iChannel,nrows,finalcol,row,channel);
 777             if((finalcol < 0) || (finalcol >= 144)) continue;
 778
 779             rocMeand->SetValue(channel,param[1]/10.0);
 780             rocRMSd->SetValue(channel,param[2]/10.0);
 781         }
 782
 783     }
 784
 785 }
 786
 787 //_______________________________________________________________________________________
 788 AliTRDCalPadStatus* AliTRDCalibPadStatus::CreateCalPadStatus()
 789 {
 790   //
 791   // Create Pad Status out of Mean and RMS values
 792   // The chamber without data are masked, this is the corrected in the preprocessor
 793   //
 794
 795   AliTRDCalPadStatus* obj = new AliTRDCalPadStatus("padstatus", "padstatus");
 796
 797   for (Int_t idet=0; idet<540; ++idet)
 798     {
 799       AliTRDCalSingleChamberStatus *calROC = obj->GetCalROC(idet);
 800
 801
 802       if ( !GetCalRocMean(idet)) {
 803         for(Int_t k = 0; k < calROC->GetNchannels(); k++){
 804           calROC->SetStatus(k,AliTRDCalPadStatus::kNotConnected);
 805         }
 806         continue;
 807       }
 808
 809
 810       //Take the stuff
 811       AliTRDCalROC *calRocMean    = new AliTRDCalROC(*( (AliTRDCalROC *) GetCalRocMean(idet)));
 812       AliTRDCalROC *calRocRMS     = new AliTRDCalROC(*( (AliTRDCalROC *) GetCalRocRMS(idet)));
 813
 814       //Take the stuff second chance
 815       AliTRDCalROC *calRocMeand    = new AliTRDCalROC(*( (AliTRDCalROC *) GetCalRocMeand(idet)));
 816       AliTRDCalROC *calRocRMSd     = new AliTRDCalROC(*( (AliTRDCalROC *) GetCalRocRMSd(idet)));
 817
 818       calRocRMS->Unfold();
 819       calRocRMSd->Unfold();
 820
 821
 822       //Range
 823       Int_t row      = calROC->GetNrows();
 824       Int_t col      = calROC->GetNcols();
 825
 826       Double_t rmsmean       = calRocMean->GetRMS()*10.0;
 827       Double_t meanmean      = calRocMean->GetMean()*10.0;
 828       Double_t meansquares   = calRocRMS->GetMean();
 829
 830
 831       for(Int_t irow = 0; irow < row; irow++){
 832
 833         // for bridged pads
 834         Float_t meanprevious = 0.0;
 835         Float_t rmsprevious  = 0.0;
 836         Float_t mean         = 0.0;
 837         Float_t rms          = 0.0;
 838
 839         for(Int_t icol = 0; icol < col; icol++){
 840
 841           mean     = calRocMean->GetValue(icol,irow)*10.0;
 842           rms      = calRocRMS->GetValue(icol,irow);
 843
 844           if(icol > 0) {
 845             meanprevious     = calRocMean->GetValue((icol -1),irow)*10.0;
 846             rmsprevious      = calRocRMS->GetValue((icol - 1),irow);
 847           }
 848
 849           Bool_t pb = kFALSE;
 850           // masked if two noisy
 851           if((rms <= 0.0001) || (TMath::Abs(mean-meanmean)>(5*rmsmean)) || (TMath::Abs(rms)>(5.0*TMath::Abs(meansquares)))) {
 852
 853             pb = kTRUE;
 854             // look at second chance
 855             Float_t meand     = calRocMeand->GetValue(icol,irow)*10.0;
 856             Float_t rmsd      = calRocRMSd->GetValue(icol,irow);
 857
 858             if((rmsd <= 0.0001) || (TMath::Abs(meand-meanmean)>(5*rmsmean)) || (TMath::Abs(rmsd)>(5.0*TMath::Abs(meansquares)))) {
 859               if((rmsd <= 0.0001) && (rms <= 0.0001)) {
 860                 calROC->SetStatus(icol,irow,AliTRDCalPadStatus::kNotConnected);
 861               }
 862               else {
 863                 calROC->SetStatus(icol, irow, AliTRDCalPadStatus::kMasked);
 864               }
 865             }
 866             else {
 867               calROC->SetStatus(icol, irow, AliTRDCalPadStatus::kReadSecond);
 868             }
 869           }
 870
 871
 872           // bridge if previous pad found something
 873           if(!pb) {
 874             if((meanprevious == mean) && (rmsprevious == rms) && (mean > 0.0001)) {
 875               //printf("mean previous %f, mean %f, rms %f, rmsprevious %f, col %d\n",meanprevious,mean,rms,rmsprevious,icol);
 876               calROC->SetStatus(icol -1 ,irow, AliTRDCalPadStatus::kPadBridgedRight);
 877               calROC->SetStatus(icol ,irow, AliTRDCalPadStatus::kPadBridgedLeft);
 878             }
 879           }
 880
 881         }
 882       }
 883
 884       delete calRocMean;
 885       delete calRocRMS;
 886       delete calRocMeand;
 887       delete calRocRMSd;
 888
 889
 890     }
 891
 892   return obj;
 893
 894 }
 895 //_______________________________________________________________________________________
 896 AliTRDCalPad* AliTRDCalibPadStatus::CreateCalPad()
 897 {
 898   //
 899   // Create Pad Noise out of RMS values
 900   //
 901
 902   AliTRDCalPad* obj = new AliTRDCalPad("PadNoise", "PadNoise");
 903
 904
 905   for (Int_t det=0; det<AliTRDgeometry::kNdet; ++det)  {
 906
 907     AliTRDCalROC *calROC22 = obj->GetCalROC(det);
 908
 909     AliTRDCalROC *calRocRMS     = ((AliTRDCalROC *)GetCalRocRMS(det,kTRUE));
 910
 911     for(Int_t k = 0; k < calROC22->GetNchannels(); k++){
 912       calROC22->SetValue(k,calRocRMS->GetValue(k));
 913     }
 914
 915   }
 916
 917   return obj;
 918
 919 }
 920
 921 //_______________________________________________________________________________________
 922 AliTRDCalDet* AliTRDCalibPadStatus::CreateCalDet() const
 923 {
 924   //
 925   // Create Det Noise correction factor
 926   //
 927
 928   AliTRDCalDet* obj = new AliTRDCalDet("DetNoise", "DetNoise (correction factor)");
 929
 930   for(Int_t l = 0; l < 540; l++){
 931     obj->SetValue(l,10.0);
 932   }
 933
 934   return obj;
 935
 936 }
 937
 938 //_____________________________________________________________________
 939 void AliTRDCalibPadStatus::DumpToFile(const Char_t *filename, const Char_t *dir, Bool_t append) /*FOLD00*/
 940 {
 941     //
 942     //  Write class to file
 943     //
 944
 945     TString sDir(dir);
 946     TString option;
 947
 948     if ( append )
 949         option = "update";
 950     else
 951         option = "recreate";
 952
 953     TDirectory *backup = gDirectory;
 954     TFile f(filename,option.Data());
 955     f.cd();
 956     if ( !sDir.IsNull() ){
 957         f.mkdir(sDir.Data());
 958         f.cd(sDir);
 959     }
 960     this->Write();
 961     f.Close();
 962
 963     if ( backup ) backup->cd();
 964 }
 965
 966 //_____________________________________________________________________
 967 void AliTRDCalibPadStatus::SetCalRocMean(AliTRDCalROC *mean, Int_t det) /*FOLD00*/
 968 {
 969     //
 970     //  Put the AliTRDCalROC in the array fCalRocArrayMean
 971     //
 972
 973
 974   AliTRDCalROC *rocMean = GetCalRocMean(det,kTRUE);
 975
 976   Int_t nChannels = rocMean->GetNchannels();
 977
 978   for (Int_t iChannel=0; iChannel<nChannels; iChannel++){
 979
 980     rocMean->SetValue(iChannel,mean->GetValue(iChannel));
 981
 982   }
 983
 984 }
 985
 986 //_____________________________________________________________________
 987 void AliTRDCalibPadStatus::SetCalRocRMS(AliTRDCalROC *rms, Int_t det) /*FOLD00*/
 988 {
 989     //
 990     //  Put the AliTRDCalROC in the array fCalRocArrayRMS
 991     //
 992
 993
 994   AliTRDCalROC *rocRms = GetCalRocRMS(det,kTRUE);
 995
 996   Int_t nChannels = rocRms->GetNchannels();
 997
 998   for (Int_t iChannel=0; iChannel<nChannels; iChannel++){
 999
1000     rocRms->SetValue(iChannel,rms->GetValue(iChannel));
1001
1002   }
1003
1004 }
1005 //_____________________________________________________________________
1006 void AliTRDCalibPadStatus::SetCalRocMeand(AliTRDCalROC *mean, Int_t det) /*FOLD00*/
1007 {
1008     //
1009     //  Put the AliTRDCalROC in the array fCalRocArrayMean
1010     //
1011
1012
1013   AliTRDCalROC *rocMean = GetCalRocMeand(det,kTRUE);
1014
1015   Int_t nChannels = rocMean->GetNchannels();
1016
1017   for (Int_t iChannel=0; iChannel<nChannels; iChannel++){
1018
1019     rocMean->SetValue(iChannel,mean->GetValue(iChannel));
1020
1021   }
1022
1023 }
1024
1025 //_____________________________________________________________________
1026 void AliTRDCalibPadStatus::SetCalRocRMSd(AliTRDCalROC *rms, Int_t det) /*FOLD00*/
1027 {
1028     //
1029     //  Put the AliTRDCalROC in the array fCalRocArrayRMS
1030     //
1031
1032
1033   AliTRDCalROC *rocRms = GetCalRocRMSd(det,kTRUE);
1034
1035   Int_t nChannels = rocRms->GetNchannels();
1036
1037   for (Int_t iChannel=0; iChannel<nChannels; iChannel++){
1038
1039     rocRms->SetValue(iChannel,rms->GetValue(iChannel));
1040
1041   }
1042
1043 }
1044 //_____________________________________________________________________________
1045 Int_t AliTRDCalibPadStatus::GetLayer(Int_t d) const
1046 {
1047   //
1048   // Reconstruct the layer number from the detector number
1049   //
1050
1051   return ((Int_t) (d % 6));
1052
1053 }
1054
1055 //_____________________________________________________________________________
1056 Int_t AliTRDCalibPadStatus::GetStack(Int_t d) const
1057 {
1058   //
1059   // Reconstruct the chamber number from the detector number
1060   //
1061
1062   return ((Int_t) (d % 30) / 6);
1063
1064 }
1065
1066 //_____________________________________________________________________________
1067 Int_t AliTRDCalibPadStatus::GetSector(Int_t d) const
1068 {
1069   //
1070   // Reconstruct the sector number from the detector number
1071   //
1072
1073   return ((Int_t) (d / 30));
1074
1075 }
1076
1077