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 **************************************************************************/
27 #include <TDirectory.h>
30 #include "AliRawReader.h"
31 #include "AliRawReaderRoot.h"
32 #include "AliRawReaderDate.h"
33 #include "AliTPCRawStream.h"
34 #include "AliTPCCalROC.h"
35 #include "AliTPCROC.h"
36 #include "AliMathBase.h"
37 #include "TTreeStream.h"
43 #include "AliTPCCalibPedestal.h"
46 ///////////////////////////////////////////////////////////////////////////////////////
47 // Implementation of the TPC pedestal and noise calibration
49 // Origin: Jens Wiechula, Marian Ivanov J.Wiechula@gsi.de, Marian.Ivanov@cern.ch
52 // *************************************************************************************
53 // * Class Description *
54 // *************************************************************************************
58 // Raw pedestal data is processed by calling one of the ProcessEvent(...) functions
59 // (see below). These in the end call the Update(...) function, where the data is filled
62 // For each ROC one TH2F histo (ROC channel vs. ADC channel) is created when
63 // it is filled for the first time (GetHistoPedestal(ROC,kTRUE)). All histos are stored in the
64 // TObjArray fHistoPedestalArray.
66 // For a fast filling of the histogram the corresponding bin number of the channel and ADC channel
67 // is computed by hand and the histogram array is accessed directly via its pointer.
68 // ATTENTION: Doing so the the entry counter of the histogram is not increased
69 // this means that e.g. the colz draw option gives an empty plot unless
70 // calling 'histo->SetEntries(1)' before drawing.
72 // After accumulating the desired statistics the Analyse() function has to be called.
73 // Whithin this function the pedestal and noise values are calculated for each pad, using
74 // the fast gaus fit function AliMathBase::FitGaus(...), and the calibration
75 // storage classes (AliTPCCalROC) are filled for each ROC.
76 // The calibration information is stored in the TObjArrays fCalRocArrayPedestal and fCalRocArrayRMS;
80 // User interface for filling data:
81 // --------------------------------
83 // To Fill information one of the following functions can be used:
85 // Bool_t ProcessEvent(eventHeaderStruct *event);
86 // - process Date event
87 // - use AliTPCRawReaderDate and call ProcessEvent(AliRawReader *rawReader)
89 // Bool_t ProcessEvent(AliRawReader *rawReader);
90 // - process AliRawReader event
91 // - use AliTPCRawStream to loop over data and call ProcessEvent(AliTPCRawStream *rawStream)
93 // Bool_t ProcessEvent(AliTPCRawStream *rawStream);
94 // - process event from AliTPCRawStream
95 // - call Update function for signal filling
97 // Int_t Update(const Int_t isector, const Int_t iRow, const Int_t
98 // iPad, const Int_t iTimeBin, const Float_t signal);
99 // - directly fill signal information (sector, row, pad, time bin, pad)
100 // to the reference histograms
102 // It is also possible to merge two independently taken calibrations using the function
104 // void Merge(AliTPCCalibPedestal *ped)
105 // - copy histograms in 'ped' if the do not exist in this instance
106 // - Add histograms in 'ped' to the histograms in this instance if the allready exist
107 // - After merging call Analyse again!
111 // -- example: filling data using root raw data:
112 // void fillPedestal(Char_t *filename)
114 // rawReader = new AliRawReaderRoot(fileName);
115 // if ( !rawReader ) return;
116 // AliTPCCalibPedestal *calib = new AliTPCCalibPedestal;
117 // while (rawReader->NextEvent()){
118 // calib->ProcessEvent(rawReader);
121 // calib->DumpToFile("PedestalData.root");
127 // What kind of information is stored and how to retrieve them:
128 // ------------------------------------------------------------
130 // - Accessing the 'Reference Histograms' (pedestal distribution histograms):
132 // TH2F *GetHistoPedestal(Int_t sector);
134 // - Accessing the calibration storage objects:
136 // AliTPCCalROC *GetCalRocPedestal(Int_t sector); - for the pedestal values
137 // AliTPCCalROC *GetCalRocNoise(Int_t sector); - for the Noise values
139 // example for visualisation:
140 // if the file "PedestalData.root" was created using the above example one could do the following:
142 // TFile filePedestal("PedestalData.root")
143 // AliTPCCalibPedestal *ped = (AliTPCCalibPedestal*)filePedestal->Get("AliTPCCalibPedestal");
144 // ped->GetCalRocPedestal(0)->Draw("colz");
145 // ped->GetCalRocRMS(0)->Draw("colz");
147 // or use the AliTPCCalPad functionality:
148 // AliTPCCalPad padPedestal(ped->GetCalPadPedestal());
149 // AliTPCCalPad padNoise(ped->GetCalPadRMS());
150 // padPedestal->MakeHisto2D()->Draw("colz"); //Draw A-Side Pedestal Information
151 // padNoise->MakeHisto2D()->Draw("colz"); //Draw A-Side Noise Information
154 example: fill pedestal with gausschen noise
155 AliTPCCalibPedestal ped;
159 TCanvas* c1 = new TCanvas;
162 ped.GetHistoPedestal(0)->SetEntries(1); //needed in order for colz to work, reason: fast filling does not increase the entries counter
163 ped.GetHistoPedestal(0)->Draw("colz");
165 ped.GetHistoPedestal(36)->SetEntries(1); //needed in order for colz to work, reason: fast filling does not increase the entries counter
166 ped.GetHistoPedestal(36)->Draw("colz");
167 TCanvas* c2 = new TCanvas;
170 ped.GetCalRocPedestal(0)->Draw("colz");
172 ped.GetCalRocRMS(0)->Draw("colz");
174 ped.GetCalRocPedestal(36)->Draw("colz");
176 ped.GetCalRocRMS(36)->Draw("colz");
183 ClassImp(AliTPCCalibPedestal)
185 AliTPCCalibPedestal::AliTPCCalibPedestal() : /*FOLD00*/
191 fOldRCUformat(kTRUE),
192 fROC(AliTPCROC::Instance()),
193 fCalRocArrayPedestal(72),
195 fHistoPedestalArray(72)
198 // default constructor
201 //_____________________________________________________________________
202 AliTPCCalibPedestal::AliTPCCalibPedestal(const AliTPCCalibPedestal &ped) : /*FOLD00*/
204 fFirstTimeBin(ped.GetFirstTimeBin()),
205 fLastTimeBin(ped.GetLastTimeBin()),
206 fAdcMin(ped.GetAdcMin()),
207 fAdcMax(ped.GetAdcMax()),
208 fOldRCUformat(kTRUE),
209 fROC(AliTPCROC::Instance()),
210 fCalRocArrayPedestal(72),
212 fHistoPedestalArray(72)
217 for (Int_t iSec = 0; iSec < 72; ++iSec){
218 const AliTPCCalROC *calPed = (AliTPCCalROC*)ped.fCalRocArrayPedestal.UncheckedAt(iSec);
219 const AliTPCCalROC *calRMS = (AliTPCCalROC*)ped.fCalRocArrayRMS.UncheckedAt(iSec);
220 const TH2F *hPed = (TH2F*)ped.fHistoPedestalArray.UncheckedAt(iSec);
222 if ( calPed != 0x0 ) fCalRocArrayPedestal.AddAt(new AliTPCCalROC(*calPed), iSec);
223 if ( calRMS != 0x0 ) fCalRocArrayRMS.AddAt(new AliTPCCalROC(*calRMS), iSec);
226 TH2F *hNew = new TH2F(*hPed);
227 hNew->SetDirectory(0);
228 fHistoPedestalArray.AddAt(hNew,iSec);
232 //_____________________________________________________________________
233 AliTPCCalibPedestal& AliTPCCalibPedestal::operator = (const AliTPCCalibPedestal &source)
236 // assignment operator
238 if (&source == this) return *this;
239 new (this) AliTPCCalibPedestal(source);
243 //_____________________________________________________________________
244 AliTPCCalibPedestal::~AliTPCCalibPedestal() /*FOLD00*/
250 fCalRocArrayPedestal.Delete();
251 fCalRocArrayRMS.Delete();
252 fHistoPedestalArray.Delete();
254 //_____________________________________________________________________
255 Int_t AliTPCCalibPedestal::Update(const Int_t icsector, /*FOLD00*/
258 const Int_t icTimeBin,
259 const Float_t csignal)
262 // Signal filling methode
265 //return if we are out of the specified time bin or adc range
266 if ( (icTimeBin>fLastTimeBin) || (icTimeBin<fFirstTimeBin) ) return 0;
267 if ( ((Int_t)csignal>fAdcMax) || ((Int_t)csignal<fAdcMin) ) return 0;
269 Int_t iChannel = fROC->GetRowIndexes(icsector)[icRow]+icPad; // global pad position in sector
271 // fast filling methode.
272 // Attention: the entry counter of the histogram is not increased
273 // this means that e.g. the colz draw option gives an empty plot
274 Int_t bin = (iChannel+1)*(fAdcMax-fAdcMin+2)+((Int_t)csignal-fAdcMin+1);
276 GetHistoPedestal(icsector,kTRUE)->GetArray()[bin]++;
280 //_____________________________________________________________________
281 Bool_t AliTPCCalibPedestal::ProcessEvent(AliTPCRawStream *rawStream)
284 // Event Processing loop - AliTPCRawStream
287 rawStream->SetOldRCUFormat(fOldRCUformat);
289 Bool_t withInput = kFALSE;
291 while (rawStream->Next()) {
293 Int_t isector = rawStream->GetSector(); // current sector
294 Int_t iRow = rawStream->GetRow(); // current row
295 Int_t iPad = rawStream->GetPad(); // current pad
296 Int_t iTimeBin = rawStream->GetTime(); // current time bin
297 Float_t signal = rawStream->GetSignal(); // current ADC signal
299 Update(isector,iRow,iPad,iTimeBin,signal);
305 //_____________________________________________________________________
306 Bool_t AliTPCCalibPedestal::ProcessEvent(AliRawReader *rawReader)
309 // Event processing loop - AliRawReader
313 AliTPCRawStream rawStream(rawReader);
315 rawReader->Select("TPC");
317 return ProcessEvent(&rawStream);
319 //_____________________________________________________________________
320 Bool_t AliTPCCalibPedestal::ProcessEvent(eventHeaderStruct *event)
323 // process date event
325 AliRawReader *rawReader = new AliRawReaderDate((void*)event);
326 Bool_t result=ProcessEvent(rawReader);
330 //_____________________________________________________________________
331 Bool_t AliTPCCalibPedestal::TestEvent() /*FOLD00*/
335 // fill one oroc and one iroc with random gaus
340 for (UInt_t iSec=0; iSec<72; ++iSec){
341 if (iSec%36>0) continue;
342 for (UInt_t iRow=0; iRow < fROC->GetNRows(iSec); ++iRow){
343 for (UInt_t iPad=0; iPad < fROC->GetNPads(iSec,iRow); ++iPad){
344 for (UInt_t iTimeBin=0; iTimeBin<1024; ++iTimeBin){
345 Float_t signal=(Int_t)(iRow+3+gRandom->Gaus(0,.7));
346 if ( signal>0 )Update(iSec,iRow,iPad,iTimeBin,signal);
353 //_____________________________________________________________________
354 TH2F* AliTPCCalibPedestal::GetHisto(Int_t sector, TObjArray *arr, /*FOLD00*/
355 Int_t nbinsY, Float_t ymin, Float_t ymax,
356 Char_t *type, Bool_t force)
359 // return pointer to Q histogram
360 // if force is true create a new histogram if it doesn't exist allready
362 if ( !force || arr->UncheckedAt(sector) )
363 return (TH2F*)arr->UncheckedAt(sector);
365 // if we are forced and histogram doesn't yes exist create it
366 Char_t name[255], title[255];
368 sprintf(name,"hCalib%s%.2d",type,sector);
369 sprintf(title,"%s calibration histogram sector %.2d;ADC channel;Channel (pad)",type,sector);
371 // new histogram with Q calib information. One value for each pad!
372 TH2F* hist = new TH2F(name,title,
374 fROC->GetNChannels(sector),0,fROC->GetNChannels(sector)
376 hist->SetDirectory(0);
377 arr->AddAt(hist,sector);
380 //_____________________________________________________________________
381 TH2F* AliTPCCalibPedestal::GetHistoPedestal(Int_t sector, Bool_t force) /*FOLD00*/
384 // return pointer to T0 histogram
385 // if force is true create a new histogram if it doesn't exist allready
387 TObjArray *arr = &fHistoPedestalArray;
388 return GetHisto(sector, arr, fAdcMax-fAdcMin, fAdcMin, fAdcMax, "Pedestal", force);
390 //_____________________________________________________________________
391 AliTPCCalROC* AliTPCCalibPedestal::GetCalRoc(Int_t sector, TObjArray* arr, Bool_t force) /*FOLD00*/
394 // return pointer to ROC Calibration
395 // if force is true create a new histogram if it doesn't exist allready
397 if ( !force || arr->UncheckedAt(sector) )
398 return (AliTPCCalROC*)arr->UncheckedAt(sector);
400 // if we are forced and the histogram doesn't yet exist create it
402 // new AliTPCCalROC for T0 information. One value for each pad!
403 AliTPCCalROC *croc = new AliTPCCalROC(sector);
404 arr->AddAt(croc,sector);
407 //_____________________________________________________________________
408 AliTPCCalROC* AliTPCCalibPedestal::GetCalRocPedestal(Int_t sector, Bool_t force) /*FOLD00*/
411 // return pointer to Carge ROC Calibration
412 // if force is true create a new histogram if it doesn't exist allready
414 TObjArray *arr = &fCalRocArrayPedestal;
415 return GetCalRoc(sector, arr, force);
417 //_____________________________________________________________________
418 AliTPCCalROC* AliTPCCalibPedestal::GetCalRocRMS(Int_t sector, Bool_t force) /*FOLD00*/
421 // return pointer to signal width ROC Calibration
422 // if force is true create a new histogram if it doesn't exist allready
424 TObjArray *arr = &fCalRocArrayRMS;
425 return GetCalRoc(sector, arr, force);
427 //_____________________________________________________________________
428 void AliTPCCalibPedestal::Merge(AliTPCCalibPedestal *ped)
431 // Merge reference histograms of sig to the current AliTPCCalibSignal
435 for (Int_t iSec=0; iSec<72; ++iSec){
436 TH2F *hRefPedMerge = ped->GetHistoPedestal(iSec);
440 TDirectory *dir = hRefPedMerge->GetDirectory(); hRefPedMerge->SetDirectory(0);
441 TH2F *hRefPed = GetHistoPedestal(iSec);
442 if ( hRefPed ) hRefPed->Add(hRefPedMerge);
444 TH2F *hist = new TH2F(*hRefPedMerge);
445 hist->SetDirectory(0);
446 fHistoPedestalArray.AddAt(hist, iSec);
448 hRefPedMerge->SetDirectory(dir);
452 //_____________________________________________________________________
453 void AliTPCCalibPedestal::Analyse() /*FOLD00*/
456 // Calculate calibration constants
459 Int_t nbinsAdc = fAdcMax-fAdcMin;
467 for (Int_t iSec=0; iSec<72; ++iSec){
468 TH2F *hP = GetHistoPedestal(iSec);
471 AliTPCCalROC *rocPedestal = GetCalRocPedestal(iSec,kTRUE);
472 AliTPCCalROC *rocRMS = GetCalRocRMS(iSec,kTRUE);
474 array_hP = hP->GetArray();
475 UInt_t nChannels = fROC->GetNChannels(iSec);
477 for (UInt_t iChannel=0; iChannel<nChannels; ++iChannel){
478 Int_t offset = (nbinsAdc+2)*(iChannel+1)+1;
479 Double_t ret = AliMathBase::FitGaus(array_hP+offset,nbinsAdc,fAdcMin,fAdcMax,¶m,&dummy);
480 // if the fitting failed set noise and pedestal to 0
485 rocPedestal->SetValue(iChannel,param[1]);
486 rocRMS->SetValue(iChannel,param[2]);
490 //_____________________________________________________________________
491 void AliTPCCalibPedestal::DumpToFile(const Char_t *filename, const Char_t *dir, Bool_t append) /*FOLD00*/
494 // Write class to file
505 TDirectory *backup = gDirectory;
506 TFile f(filename,option.Data());
508 if ( !sDir.IsNull() ){
509 f.mkdir(sDir.Data());
515 if ( backup ) backup->cd();