TPC/AliTPCdataQA.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
  17 /* $Id$ */
  18
  19
  20 // stl includes
  21 #include <iostream>
  22
  23 using namespace std;
  24
  25 //Root includes
  26 #include <TH1F.h>
  27 #include <TH2F.h>
  28 #include <TString.h>
  29 #include <TMath.h>
  30 #include <TF1.h>
  31 #include <TRandom.h>
  32 #include <TDirectory.h>
  33 #include <TFile.h>
  34 //AliRoot includes
  35 #include "AliRawReader.h"
  36 #include "AliRawReaderRoot.h"
  37 #include "AliRawReaderDate.h"
  38 #include "AliTPCRawStream.h"
  39 #include "AliTPCCalROC.h"
  40 #include "AliTPCROC.h"
  41 #include "AliMathBase.h"
  42 #include "TTreeStream.h"
  43 #include "AliTPCRawStreamFast.h"
  44
  45 //date
  46 #include "event.h"
  47 #include "AliTPCCalPad.h"
  48
  49 //header file
  50 #include "AliTPCdataQA.h"
  51
  52 ClassImp(AliTPCdataQA)
  53
  54 AliTPCdataQA::AliTPCdataQA() : /*FOLD00*/
  55   TObject(),
  56   fFirstTimeBin(60),
  57   fLastTimeBin(1000),
  58   fAdcMin(1),
  59   fAdcMax(100),
  60   fOldRCUformat(kTRUE),
  61   fROC(AliTPCROC::Instance()),
  62   fMapping(NULL),
  63   fPedestal(0),
  64   fNoise(0),
  65   fMaxCharge(0),
  66   fMeanCharge(0),
  67   fNoThreshold(0),
  68   fOverThreshold0(0),
  69   fOverThreshold5(0),
  70   fOverThreshold10(0),
  71   fOverThreshold20(0),
  72   fOverThreshold30(0),
  73   fEventCounter(0)
  74 {
  75   //
  76   // default constructor
  77   //
  78
  79   fSectorLast  = -1;
  80   fRowLast     =  0;
  81   fPadLast     =  0;
  82   fTimeBinLast =  0;
  83   fSignalLast  =  0;
  84   fNAboveThreshold = 0;
  85 }
  86
  87
  88 //_____________________________________________________________________
  89 AliTPCdataQA::AliTPCdataQA(const AliTPCdataQA &ped) : /*FOLD00*/
  90   TObject(ped),
  91   fFirstTimeBin(ped.GetFirstTimeBin()),
  92   fLastTimeBin(ped.GetLastTimeBin()),
  93   fAdcMin(ped.GetAdcMin()),
  94   fAdcMax(ped.GetAdcMax()),
  95   fOldRCUformat(ped.fOldRCUformat),
  96   fROC(AliTPCROC::Instance()),
  97   fMapping(NULL)
  98 {
  99   //
 100   // copy constructor
 101   //
 102
 103 }
 104
 105
 106 //_____________________________________________________________________
 107 AliTPCdataQA& AliTPCdataQA::operator = (const  AliTPCdataQA &source)
 108 {
 109   //
 110   // assignment operator
 111   //
 112   if (&source == this) return *this;
 113   new (this) AliTPCdataQA(source);
 114
 115   return *this;
 116 }
 117
 118
 119 //_____________________________________________________________________
 120 AliTPCdataQA::~AliTPCdataQA() /*FOLD00*/
 121 {
 122   //
 123   // destructor
 124   //
 125
 126   // do not delete fMapping, because we do not own it.
 127
 128 }
 129
 130
 131
 132
 133 //_____________________________________________________________________
 134 Bool_t AliTPCdataQA::ProcessEventFast(AliTPCRawStreamFast *rawStreamFast)
 135 {
 136   //
 137   // Event Processing loop - AliTPCRawStream
 138   //
 139   Bool_t withInput = kFALSE;
 140
 141   while ( rawStreamFast->NextDDL() ){
 142       while ( rawStreamFast->NextChannel() ){
 143           Int_t isector  = rawStreamFast->GetSector();                       //  current sector
 144           Int_t iRow     = rawStreamFast->GetRow();                          //  current row
 145           Int_t iPad     = rawStreamFast->GetPad();                          //  current pad
 146
 147           while ( rawStreamFast->NextBunch() ){
 148             Int_t startTbin = (Int_t)rawStreamFast->GetStartTimeBin();
 149             Int_t endTbin = (Int_t)rawStreamFast->GetEndTimeBin();
 150
 151               for (Int_t iTimeBin = startTbin; iTimeBin < endTbin; iTimeBin++){
 152                   Float_t signal=(Float_t)rawStreamFast->GetSignals()[iTimeBin-startTbin];
 153                   Update(isector,iRow,iPad,iTimeBin+1,signal);
 154                   withInput = kTRUE;
 155               }
 156           }
 157       }
 158   }
 159
 160   return withInput;
 161 }
 162 //_____________________________________________________________________
 163 Bool_t AliTPCdataQA::ProcessEventFast(AliRawReader *rawReader)
 164 {
 165   //
 166   //  Event processing loop - AliRawReader
 167   //
 168   fSectorLast  = -1;
 169   AliTPCRawStreamFast *rawStreamFast = new AliTPCRawStreamFast(rawReader, (AliAltroMapping**)fMapping);
 170   Bool_t res=ProcessEventFast(rawStreamFast);
 171   if(res)
 172     fEventCounter++; // only increment event counter if there is TPC data
 173                      // otherwise Analyse (called in QA) fails
 174
 175   delete rawStreamFast;
 176   return res;
 177 }
 178
 179 //_____________________________________________________________________
 180 Bool_t AliTPCdataQA::ProcessEvent(AliTPCRawStream *rawStream)
 181 {
 182   //
 183   // Event Processing loop - AliTPCRawStream
 184   //
 185
 186   rawStream->SetOldRCUFormat(fOldRCUformat);
 187
 188   Bool_t withInput = kFALSE;
 189
 190   while (rawStream->Next()) {
 191
 192     Int_t iSector  = rawStream->GetSector();      //  current ROC
 193     Int_t iRow     = rawStream->GetRow();         //  current row
 194     Int_t iPad     = rawStream->GetPad();         //  current pad
 195     Int_t iTimeBin = rawStream->GetTime();        //  current time bin
 196     Float_t signal = rawStream->GetSignal();      //  current ADC signal
 197
 198     Update(iSector,iRow,iPad,iTimeBin,signal);
 199     withInput = kTRUE;
 200   }
 201
 202   return withInput;
 203 }
 204
 205
 206 //_____________________________________________________________________
 207 Bool_t AliTPCdataQA::ProcessEvent(AliRawReader *rawReader)
 208 {
 209   //
 210   //  Event processing loop - AliRawReader
 211   //
 212
 213   // if fMapping is NULL the rawstream will crate its own mapping
 214   fSectorLast  = -1;
 215   AliTPCRawStream rawStream(rawReader, (AliAltroMapping**)fMapping);
 216   rawReader->Select("TPC");
 217   Bool_t res =  ProcessEvent(&rawStream);
 218
 219   if(res)
 220     fEventCounter++; // only increment event counter if there is TPC data
 221                      // otherwise Analyse (called in QA) fails
 222   return res;
 223 }
 224
 225
 226 //_____________________________________________________________________
 227 Bool_t AliTPCdataQA::ProcessEvent(eventHeaderStruct *event)
 228 {
 229   //
 230   //  process date event
 231   //
 232
 233   AliRawReader *rawReader = new AliRawReaderDate((void*)event);
 234   Bool_t result=ProcessEvent(rawReader);
 235   delete rawReader;
 236   return result;
 237 }
 238
 239
 240
 241 //_____________________________________________________________________
 242 void AliTPCdataQA::DumpToFile(const Char_t *filename, const Char_t *dir, Bool_t append) /*FOLD00*/
 243 {
 244   //
 245   //  Write class to file
 246   //
 247
 248   TString sDir(dir);
 249   TString option;
 250
 251   if ( append )
 252     option = "update";
 253   else
 254     option = "recreate";
 255
 256   TDirectory *backup = gDirectory;
 257   TFile f(filename,option.Data());
 258   f.cd();
 259   if ( !sDir.IsNull() ){
 260     f.mkdir(sDir.Data());
 261     f.cd(sDir);
 262   }
 263   this->Write();
 264   f.Close();
 265
 266   if ( backup ) backup->cd();
 267 }
 268
 269
 270 //_____________________________________________________________________
 271 Int_t AliTPCdataQA::Update(const Int_t icsector, /*FOLD00*/
 272                                   const Int_t icRow,
 273                                   const Int_t icPad,
 274                                   const Int_t icTimeBin,
 275                                   const Float_t csignal)
 276 {
 277   //
 278   // Signal filling method
 279   //
 280   if (icTimeBin<fFirstTimeBin) return 0;
 281   if (icTimeBin>fLastTimeBin) return 0;
 282
 283   if (!fMaxCharge) fMaxCharge = new AliTPCCalPad("MaxCharge","MaxCharge");
 284   if (!fMeanCharge) fMeanCharge = new AliTPCCalPad("MeanCharge","MeanCharge");
 285   if (!fNoThreshold) fNoThreshold = new AliTPCCalPad("NoThreshold","NoThreshold");
 286   if (!fOverThreshold0) fOverThreshold0 = new AliTPCCalPad("OverThreshold0","OverThreshold0");
 287   if (!fOverThreshold5) fOverThreshold5 = new AliTPCCalPad("OverThreshold5","OverThreshold5");
 288   if (!fOverThreshold10) fOverThreshold10 = new AliTPCCalPad("OverThreshold10","OverThreshold10");
 289   if (!fOverThreshold20) fOverThreshold20 = new AliTPCCalPad("OverThreshold20","OverThreshold20");
 290   if (!fOverThreshold30) fOverThreshold30 = new AliTPCCalPad("OverThreshold30","OverThreshold30");
 291   //
 292
 293   Int_t signal = Int_t(csignal);
 294
 295   // if pedestal calibrations are loaded subtract pedestals
 296   if(fPedestal) {
 297
 298     Int_t pedestal = Int_t(fPedestal->GetCalROC(icsector)->GetValue(icRow, icPad));
 299     if(pedestal<10 || pedestal>90)
 300       return 0;
 301     signal -= pedestal;
 302   }
 303
 304
 305   if (signal >= 0) {
 306
 307     Float_t count = fNoThreshold->GetCalROC(icsector)->GetValue(icRow, icPad);
 308     fNoThreshold->GetCalROC(icsector)->SetValue(icRow, icPad,count+1);
 309   }
 310
 311   // Require at least 3 ADC channels
 312   if (signal < 3)
 313     return 0;
 314
 315   // if noise calibrations are loaded require at least 3*sigmaNoise
 316   if(fNoise) {
 317
 318     Float_t noise = fNoise->GetCalROC(icsector)->GetValue(icRow, icPad);
 319
 320     if(signal<noise*3)
 321       return 0;
 322   }
 323   //
 324   // this signal is ok - now see if the previous signal was connected
 325   // this is a bit ugly since the standard decoder goes down in time bins
 326   // (10, 9, 8..) while the fast HLT decoder goes up in time bins (1, 2, 3..)
 327   //
 328   if(fSectorLast==icsector && fRowLast==icRow && fPadLast==icPad &&
 329      fTimeBinLast==icTimeBin+1 || fTimeBinLast==icTimeBin-1)
 330     fNAboveThreshold++;
 331   else
 332     fNAboveThreshold = 1;
 333
 334   if(fNAboveThreshold==2) {
 335
 336     //
 337     // This is the only special case, because before we did not know if we
 338     // should store the information
 339     //
 340     UpdateSignalHistograms(fSectorLast, fRowLast, fPadLast, fTimeBinLast,
 341                          fSignalLast);
 342   }
 343
 344   // keep the information for the next signal
 345   fSectorLast  = icsector;
 346   fRowLast     = icRow;
 347   fPadLast     = icPad;
 348   fTimeBinLast = icTimeBin;
 349   fSignalLast  = signal;
 350
 351   if(fNAboveThreshold==1) // we don't know if this information should be stored
 352     return 1;
 353
 354   UpdateSignalHistograms(fSectorLast, fRowLast, fPadLast, fTimeBinLast,
 355                        fSignalLast);
 356
 357   return 1;
 358 }
 359 //_____________________________________________________________________
 360 void AliTPCdataQA::UpdateSignalHistograms(const Int_t icsector, /*FOLD00*/
 361                                         const Int_t icRow,
 362                                         const Int_t icPad,
 363                                         const Int_t icTimeBin,
 364                                         const Float_t signal)
 365 {
 366   //
 367   // Signal filling method
 368   //
 369
 370   {
 371     Float_t charge = fMeanCharge->GetCalROC(icsector)->GetValue(icRow, icPad);
 372     fMeanCharge->GetCalROC(icsector)->SetValue(icRow, icPad, charge + signal);
 373   }
 374
 375   if (signal>fMaxCharge->GetCalROC(icsector)->GetValue(icRow, icPad)){
 376     fMaxCharge->GetCalROC(icsector)->SetValue(icRow, icPad,signal);
 377   }
 378
 379   if (signal>0){
 380     Float_t count = fOverThreshold0->GetCalROC(icsector)->GetValue(icRow, icPad);
 381     fOverThreshold0->GetCalROC(icsector)->SetValue(icRow, icPad,count+1);
 382   };
 383   //
 384   if (signal>5){
 385     Float_t count = fOverThreshold5->GetCalROC(icsector)->GetValue(icRow, icPad);
 386     fOverThreshold5->GetCalROC(icsector)->SetValue(icRow, icPad,count+1);
 387   };
 388   if (signal>10){
 389     Float_t count = fOverThreshold10->GetCalROC(icsector)->GetValue(icRow, icPad);
 390     fOverThreshold10->GetCalROC(icsector)->SetValue(icRow, icPad,count+1);
 391   };
 392   if (signal>20){
 393     Float_t count = fOverThreshold20->GetCalROC(icsector)->GetValue(icRow, icPad);
 394     fOverThreshold20->GetCalROC(icsector)->SetValue(icRow, icPad,count+1);
 395   };
 396   if (signal>30){
 397     Float_t count = fOverThreshold30->GetCalROC(icsector)->GetValue(icRow, icPad);
 398     fOverThreshold30->GetCalROC(icsector)->SetValue(icRow, icPad,count+1);
 399   };
 400 }
 401
 402 //_____________________________________________________________________
 403 void AliTPCdataQA::Analyse()
 404 {
 405   //
 406   //  Calculate calibration constants
 407   //
 408
 409   cout << "Analyse called" << endl;
 410
 411   if(fEventCounter==0) {
 412
 413     cout << "EventCounter == 0, Cannot analyse" << endl;
 414     return;
 415   }
 416
 417   Int_t nTimeBins = fLastTimeBin - fFirstTimeBin +1;
 418
 419   cout << "EventCounter: " << fEventCounter << endl
 420        << "TimeBins: " << nTimeBins << endl;
 421
 422   if (fMeanCharge && fNoThreshold) fMeanCharge->Divide(fNoThreshold);
 423
 424   Float_t normalization = 1.0 / Float_t(fEventCounter * nTimeBins);
 425   if (fNoThreshold)     fNoThreshold->Multiply(normalization);
 426   if (fOverThreshold0)  fOverThreshold0->Multiply(normalization);
 427   if (fOverThreshold5)  fOverThreshold5->Multiply(normalization);
 428   if (fOverThreshold10) fOverThreshold10->Multiply(normalization);
 429   if (fOverThreshold20) fOverThreshold20->Multiply(normalization);
 430   if (fOverThreshold30) fOverThreshold30->Multiply(normalization);
 431 }