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 **************************************************************************/
17 //-------------------------------------------------------
18 // Implementation of the TPC pedestal and noise calibration
20 // Origin: Jens Wiechula, Marian Ivanov J.Wiechula@gsi.de, Marian.Ivanov@cern.ch
23 //-------------------------------------------------------
29 example: fill pedestal with gausschen noise
30 AliTPCCalibPedestal ped;
34 TCanvas* c1 = new TCanvas;
37 ped.GetHistoPedestal(0)->SetEntries(1); //needed in order for colz to work, reason: fast filling does not increase the entries counter
38 ped.GetHistoPedestal(0)->Draw("colz");
40 ped.GetHistoPedestal(36)->SetEntries(1); //needed in order for colz to work, reason: fast filling does not increase the entries counter
41 ped.GetHistoPedestal(36)->Draw("colz");
42 TCanvas* c2 = new TCanvas;
45 ped.GetCalRocPedestal(0)->Draw("colz");
47 ped.GetCalRocRMS(0)->Draw("colz");
49 ped.GetCalRocPedestal(36)->Draw("colz");
51 ped.GetCalRocRMS(36)->Draw("colz");
57 #include <TObjArray.h>
64 #include <TDirectory.h>
67 #include "AliRawReader.h"
68 #include "AliRawReaderRoot.h"
69 #include "AliRawReaderDate.h"
70 #include "AliTPCRawStream.h"
71 #include "AliTPCCalROC.h"
72 #include "AliTPCROC.h"
73 #include "AliMathBase.h"
74 #include "TTreeStream.h"
80 #include "AliTPCCalibPedestal.h"
83 ClassImp(AliTPCCalibPedestal) /*FOLD00*/
91 AliTPCCalibPedestal::AliTPCCalibPedestal() : /*FOLD00*/
97 fROC(AliTPCROC::Instance()),
98 fCalRocArrayPedestal(72),
100 fHistoPedestalArray(72)
103 // default constructor
106 //_____________________________________________________________________
107 AliTPCCalibPedestal::AliTPCCalibPedestal(const AliTPCCalibPedestal &ped) : /*FOLD00*/
109 fFirstTimeBin(ped.GetFirstTimeBin()),
110 fLastTimeBin(ped.GetLastTimeBin()),
111 fAdcMin(ped.GetAdcMin()),
112 fAdcMax(ped.GetAdcMax()),
113 fROC(AliTPCROC::Instance()),
114 fCalRocArrayPedestal(72),
116 fHistoPedestalArray(72)
121 for (Int_t iSec = 0; iSec < 72; iSec++){
122 const AliTPCCalROC *calPed = (AliTPCCalROC*)ped.fCalRocArrayPedestal.UncheckedAt(iSec);
123 const AliTPCCalROC *calRMS = (AliTPCCalROC*)ped.fCalRocArrayRMS.UncheckedAt(iSec);
124 const TH2F *hPed = (TH2F*)ped.fHistoPedestalArray.UncheckedAt(iSec);
126 if ( calPed != 0x0 ) fCalRocArrayPedestal.AddAt(new AliTPCCalROC(*calPed), iSec);
127 if ( calRMS != 0x0 ) fCalRocArrayRMS.AddAt(new AliTPCCalROC(*calRMS), iSec);
130 TH2F *hNew = new TH2F(*hPed);
131 hNew->SetDirectory(0);
132 fHistoPedestalArray.AddAt(hNew,iSec);
136 //_____________________________________________________________________
137 AliTPCCalibPedestal& AliTPCCalibPedestal::operator = (const AliTPCCalibPedestal &source)
140 // assignment operator
142 if (&source == this) return *this;
143 new (this) AliTPCCalibPedestal(source);
147 //_____________________________________________________________________
148 AliTPCCalibPedestal::~AliTPCCalibPedestal() /*FOLD00*/
159 //_____________________________________________________________________
160 Int_t AliTPCCalibPedestal::Update(const Int_t icsector, /*FOLD00*/
163 const Int_t icTimeBin,
164 const Float_t csignal)
167 // Signal filling methode
169 if ( (icTimeBin>fLastTimeBin) || (icTimeBin<fFirstTimeBin) ) return 0;
171 Int_t iChannel = fROC->GetRowIndexes(icsector)[icRow]+icPad; // global pad position in sector
173 // fast filling methode.
174 // Attention: the entry counter of the histogram is not increased
175 // this means that e.g. the colz draw option gives an empty plot
177 if ( !(((Int_t)csignal>fAdcMax ) || ((Int_t)csignal<fAdcMin)) )
178 bin = (iChannel+1)*(fAdcMax-fAdcMin+2)+((Int_t)csignal-fAdcMin+1);
180 GetHistoPedestal(icsector,kTRUE)->GetArray()[bin]++;
184 //_____________________________________________________________________
185 Bool_t AliTPCCalibPedestal::ProcessEvent(AliTPCRawStream *rawStream)
188 // Event Processing loop - AliTPCRawStream
191 rawStream->SetOldRCUFormat(1);
193 Bool_t withInput = kFALSE;
195 while (rawStream->Next()) {
197 Int_t isector = rawStream->GetSector(); // current sector
198 Int_t iRow = rawStream->GetRow(); // current row
199 Int_t iPad = rawStream->GetPad(); // current pad
200 Int_t iTimeBin = rawStream->GetTime(); // current time bin
201 Float_t signal = rawStream->GetSignal(); // current ADC signal
203 Update(isector,iRow,iPad,iTimeBin,signal);
209 //_____________________________________________________________________
210 Bool_t AliTPCCalibPedestal::ProcessEvent(AliRawReader *rawReader)
213 // Event processing loop - AliRawReader
217 AliTPCRawStream rawStream(rawReader);
219 rawReader->Select("TPC");
221 return ProcessEvent(&rawStream);
223 //_____________________________________________________________________
224 Bool_t AliTPCCalibPedestal::ProcessEvent(eventHeaderStruct *event)
227 // process date event
229 AliRawReader *rawReader = new AliRawReaderDate((void*)event);
230 Bool_t result=ProcessEvent(rawReader);
234 //_____________________________________________________________________
235 Bool_t AliTPCCalibPedestal::TestEvent() /*FOLD00*/
239 // fill one oroc and one iroc with random gaus
244 for (UInt_t iSec=0; iSec<72; iSec++){
245 if (iSec%36>0) continue;
246 for (UInt_t iRow=0; iRow < fROC->GetNRows(iSec); iRow++){
247 for (UInt_t iPad=0; iPad < fROC->GetNPads(iSec,iRow); iPad++){
248 for (UInt_t iTimeBin=0; iTimeBin<1024; iTimeBin++){
249 Float_t signal=(Int_t)(iRow+3+gRandom->Gaus(0,.7));
250 if ( signal>0 )Update(iSec,iRow,iPad,iTimeBin,signal);
257 //_____________________________________________________________________
258 TH2F* AliTPCCalibPedestal::GetHisto(Int_t sector, TObjArray *arr, /*FOLD00*/
259 Int_t nbinsY, Float_t ymin, Float_t ymax,
260 Char_t *type, Bool_t force)
263 // return pointer to Q histogram
264 // if force is true create a new histogram if it doesn't exist allready
266 if ( !force || arr->UncheckedAt(sector) )
267 return (TH2F*)arr->UncheckedAt(sector);
269 // if we are forced and histogram doesn't yes exist create it
270 Char_t name[255], title[255];
272 sprintf(name,"hCalib%s%.2d",type,sector);
273 sprintf(title,"%s calibration histogram sector %.2d;ADC channel;Channel (pad)",type,sector);
275 // new histogram with Q calib information. One value for each pad!
276 TH2F* hist = new TH2F(name,title,
278 fROC->GetNChannels(sector),0,fROC->GetNChannels(sector)
280 hist->SetDirectory(0);
281 arr->AddAt(hist,sector);
284 //_____________________________________________________________________
285 TH2F* AliTPCCalibPedestal::GetHistoPedestal(Int_t sector, Bool_t force) /*FOLD00*/
288 // return pointer to T0 histogram
289 // if force is true create a new histogram if it doesn't exist allready
291 TObjArray *arr = &fHistoPedestalArray;
292 return GetHisto(sector, arr, fAdcMax-fAdcMin, fAdcMin, fAdcMax, "Pedestal", force);
294 //_____________________________________________________________________
295 AliTPCCalROC* AliTPCCalibPedestal::GetCalRoc(Int_t sector, TObjArray* arr, Bool_t force) /*FOLD00*/
298 // return pointer to ROC Calibration
299 // if force is true create a new histogram if it doesn't exist allready
301 if ( !force || arr->UncheckedAt(sector) )
302 return (AliTPCCalROC*)arr->UncheckedAt(sector);
304 // if we are forced and histogram doesn't yes exist create it
306 // new AliTPCCalROC for T0 information. One value for each pad!
307 AliTPCCalROC *croc = new AliTPCCalROC(sector);
309 for ( UInt_t iChannel = 0; iChannel<croc->GetNchannels(); iChannel++){
310 croc->SetValue(iChannel, 0);
312 arr->AddAt(croc,sector);
315 //_____________________________________________________________________
316 AliTPCCalROC* AliTPCCalibPedestal::GetCalRocPedestal(Int_t sector, Bool_t force) /*FOLD00*/
319 // return pointer to Carge ROC Calibration
320 // if force is true create a new histogram if it doesn't exist allready
322 TObjArray *arr = &fCalRocArrayPedestal;
323 return GetCalRoc(sector, arr, force);
325 //_____________________________________________________________________
326 AliTPCCalROC* AliTPCCalibPedestal::GetCalRocRMS(Int_t sector, Bool_t force) /*FOLD00*/
329 // return pointer to signal width ROC Calibration
330 // if force is true create a new histogram if it doesn't exist allready
332 TObjArray *arr = &fCalRocArrayRMS;
333 return GetCalRoc(sector, arr, force);
335 //_____________________________________________________________________
336 void AliTPCCalibPedestal::Analyse() /*FOLD00*/
339 // Calculate calibration constants
342 Int_t nbinsAdc = fAdcMax-fAdcMin;
350 for (Int_t iSec=0; iSec<72; iSec++){
351 TH2F *hP = GetHistoPedestal(iSec);
354 AliTPCCalROC *rocPedestal = GetCalRocPedestal(iSec,kTRUE);
355 AliTPCCalROC *rocRMS = GetCalRocRMS(iSec,kTRUE);
357 array_hP = hP->GetArray();
358 UInt_t nChannels = fROC->GetNChannels(iSec);
360 for (UInt_t iChannel=0; iChannel<nChannels; iChannel++){
361 Int_t offset = (nbinsAdc+2)*(iChannel+1)+1;
362 Double_t ret = AliMathBase::FitGaus(array_hP+offset,nbinsAdc,fAdcMin,fAdcMax,¶m,&dummy);
363 // if the fitting failed set noise and pedestal to 0
368 rocPedestal->SetValue(iChannel,param[1]);
369 rocRMS->SetValue(iChannel,param[2]);
373 //_____________________________________________________________________
374 void AliTPCCalibPedestal::DumpToFile(const Char_t *filename, const Char_t *dir, Bool_t append) /*FOLD00*/
377 // Write class to file
388 TDirectory *backup = gDirectory;
389 TFile f(filename,option.Data());
391 if ( !sDir.IsNull() ){
392 f.mkdir(sDir.Data());
398 if ( backup ) backup->cd();