]> git.uio.no Git - u/mrichter/AliRoot.git/blame - TPC/Calib/AliTPCCalibKr.cxx
*** V interface for TPCCalibTasks ***
[u/mrichter/AliRoot.git] / TPC / Calib / AliTPCCalibKr.cxx
CommitLineData
a65a7e70 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//Root includes
18#include <TH1F.h>
19#include <TH1D.h>
20#include <TH2F.h>
21#include <TH3F.h>
22#include <TString.h>
23#include <TMath.h>
24#include <TF1.h>
25#include <TRandom.h>
26#include <TDirectory.h>
27#include <TFile.h>
28#include <TAxis.h>
29//AliRoot includes
30#include "AliRawReader.h"
31#include "AliRawReaderRoot.h"
32#include "AliRawReaderDate.h"
33#include "AliTPCCalROC.h"
34#include "AliTPCCalPad.h"
35#include "AliTPCROC.h"
36#include "AliMathBase.h"
37#include "TTreeStream.h"
38
39//date
40#include "event.h"
41
42//header file
43#include "AliTPCCalibKr.h"
44
45//----------------------------------------------------------------------------
46// The AliTPCCalibKr class description (TPC Kr calibration).
47//
48//
49// The AliTPCCalibKr keeps the array of TH3F histograms (TPC_max_padraw,TPC_max_pad,TPC_ADC_cluster),
50// its data memebers and is filled by AliTPCCalibKrTask under conditions (Accept()).
51//
52// The ouput TH3F histograms are later used to determine the calibration parameters of TPC chambers.
53// These calculations are done by using AliTPCCalibKr::Analyse() function. The ouput calibration
54// parameters (details in AliTPCCalibKr::Analyse()) are stored in the calibKr.root file for each TPC pad.
55// In addition the debugCalibKr.root file with debug information is created.
56//
57
58/*
59
60// Usage example:
61//
62
63// 1. Analyse output histograms:
64TFile f("outHistFile.root");
65AliTPCCalibKr *obj = (AliTPCCalibKr*) cOutput.FindObject("AliTPCCalibKr");
66obj->SetRadius(0,0);
67obj->SetStep(1,1);
68obj->Analyse();
69
70// 2. See calibration parameters e.g.:
71TFile f("calibKr.root");
72spectrMean->GetCalROC(70)->GetValue(40,40);
73fitMean->GetCalROC(70)->GetValue(40,40);
74
75// 3. See debug information e.g.:
76TFile f("debugCalibKr.root");
77.ls;
78
79// -- Print calibKr TTree content
80calibKr->Print();
81
82// -- Draw calibKr TTree variables
83calibKr.Draw("fitMean");
84
85*/
86
87
88//
89// Author: Jacek Otwinowski (J.Otwinowski@gsi.de) and Stafan Geartner (S.Gaertner@gsi.de)
90//-----------------------------------------------------------------------------
91
92ClassImp(AliTPCCalibKr)
93
94AliTPCCalibKr::AliTPCCalibKr() :
95 TObject(),
96 fASide(kTRUE),
97 fCSide(kTRUE),
98 fHistoKrArray(72),
99 fADCOverClustSizeMin(0.0),
100 fADCOverClustSizeMax(1.0e9),
101 fMaxADCOverClustADCMin(0.0),
102 fMaxADCOverClustADCMax(1.0e9),
103 fTimeMin(0.0),
104 fTimeMax(1.0e9),
105 fClustSizeMin(0.0),
106 fClustSizeMax(1.0e9),
107 fTimebinRmsIrocMin(0.0),
108 fPadRmsIrocMin(0.0),
109 fRowRmsIrocMin(0.0),
110 fClusterPadSize1DIrocMax(200),
111 fCurveCoefficientIroc(1.0e9),
112 fTimebinRmsOrocMin(0.0),
113 fPadRmsOrocMin(0.0),
114 fRowRmsOrocMin(0.0),
115 fClusterPadSize1DOrocMax(200),
116 fCurveCoefficientOroc(1.0e9),
117 fIrocHistogramMin(100.),
118 fIrocHistogramMax(6000.),
119 fIrocHistogramNbins(200),
120 fOrocHistogramMin(100.),
121 fOrocHistogramMax(5500.),
122 fOrocHistogramNbins(200),
123 fRowRadius(0),
124 fPadRadius(0),
125 fRowStep(1),
126 fPadStep(1)
127
128{
129 //
130 // default constructor
131 //
132
133 // TObjArray with histograms
134 fHistoKrArray.SetOwner(kTRUE); // is owner of histograms
135 fHistoKrArray.Clear();
136
137 // init cuts (by Stefan)
138// SetADCOverClustSizeRange(7,200);
139// SetMaxADCOverClustADCRange(0.01,0.4);
140// SetTimeRange(200,600);
141// SetClustSizeRange(6,200);
142
143 //init cuts (by Adam)
144 SetTimebinRmsMin(1.6,0.8);
145 SetPadRmsMin(0.825,0.55);
146 SetRowRmsMin(0.1,0.1);
147 SetClusterPadSize1DMax(15,11);
148 SetCurveCoefficient(1.,2.);
149
150 //set histograms settings
151 SetIrocHistogram(200,100,6000);
152 SetOrocHistogram(200,100,5500);
153
154 // init histograms
155 //Init();
156}
157
158//_____________________________________________________________________
159AliTPCCalibKr::AliTPCCalibKr(const AliTPCCalibKr& pad) :
160 TObject(pad),
161
162 fASide(pad.fASide),
163 fCSide(pad.fCSide),
164 fHistoKrArray(72),
165 fADCOverClustSizeMin(pad.fADCOverClustSizeMin),
166 fADCOverClustSizeMax(pad.fADCOverClustSizeMax),
167 fMaxADCOverClustADCMin(pad.fMaxADCOverClustADCMin),
168 fMaxADCOverClustADCMax(pad.fMaxADCOverClustADCMax),
169 fTimeMin(pad.fTimeMin),
170 fTimeMax(pad.fTimeMax),
171 fClustSizeMin(pad.fClustSizeMin),
172 fClustSizeMax(pad.fClustSizeMax),
173 fTimebinRmsIrocMin(pad.fTimebinRmsIrocMin),
174 fPadRmsIrocMin(pad.fPadRmsIrocMin),
175 fRowRmsIrocMin(pad.fRowRmsIrocMin),
176 fClusterPadSize1DIrocMax(pad.fClusterPadSize1DIrocMax),
177 fCurveCoefficientIroc(pad.fCurveCoefficientIroc),
178 fTimebinRmsOrocMin(pad.fTimebinRmsOrocMin),
179 fPadRmsOrocMin(pad.fPadRmsOrocMin),
180 fRowRmsOrocMin(pad.fRowRmsOrocMin),
181 fClusterPadSize1DOrocMax(pad.fClusterPadSize1DOrocMax),
182 fCurveCoefficientOroc(pad.fCurveCoefficientOroc),
183 fIrocHistogramMin(pad.fIrocHistogramMin),
184 fIrocHistogramMax(pad.fIrocHistogramMax),
185 fIrocHistogramNbins(pad.fIrocHistogramNbins),
186 fOrocHistogramMin(pad.fOrocHistogramMin),
187 fOrocHistogramMax(pad.fOrocHistogramMax),
188 fOrocHistogramNbins(pad.fOrocHistogramNbins),
189 fRowRadius(pad.fRowRadius),
190 fPadRadius(pad.fPadRadius),
191 fRowStep(pad.fRowStep),
192 fPadStep(pad.fPadStep)
193
194{
195 // copy constructor
196
197 // TObjArray with histograms
198 fHistoKrArray.SetOwner(kTRUE); // is owner of histograms
199 fHistoKrArray.Clear();
200
201 for (Int_t iSec = 0; iSec < 72; ++iSec)
202 {
203 TH3F *hOld = pad.GetHistoKr(iSec);
204 if(hOld) {
205 TH3F *hNew = new TH3F( *pad.GetHistoKr(iSec) );
206 fHistoKrArray.AddAt(hNew,iSec);
207 }
208 }
209}
210
211//_____________________________________________________________________
212AliTPCCalibKr::~AliTPCCalibKr()
213{
214 //
215 // destructor
216 //
217
218 // for (Int_t iSec = 0; iSec < 72; ++iSec) {
219 // if (fHistoKrArray.At(iSec)) delete fHistoKrArray.RemoveAt(iSec);
220 // }
221 fHistoKrArray.Delete();
222}
223
224//_____________________________________________________________________
225AliTPCCalibKr& AliTPCCalibKr::operator = (const AliTPCCalibKr &source)
226{
227 // assignment operator
228
229 if (&source == this) return *this;
230 new (this) AliTPCCalibKr(source);
231
232 return *this;
233}
234
235//_____________________________________________________________________
236void AliTPCCalibKr::Init()
237{
238 //
239 // init output histograms
240 //
241
242 // add histograms to the TObjArray
243 for(Int_t i=0; i<72; ++i) {
244
245 // C - side
246 if(IsCSide(i) == kTRUE && fCSide == kTRUE) {
247 TH3F *hist = CreateHisto(i);
248 if(hist) fHistoKrArray.AddAt(hist,i);
249 }
250
251 // A - side
252 if(IsCSide(i) == kFALSE && fASide == kTRUE) {
253 TH3F *hist = CreateHisto(i);
254 if(hist) fHistoKrArray.AddAt(hist,i);
255 }
256 }
257}
258
259//_____________________________________________________________________
260Bool_t AliTPCCalibKr::Process(AliTPCclusterKr *cluster)
261{
262 //
263 // process events
264 // call event by event
265 //
266
267 if(cluster) Update(cluster);
268 else return kFALSE;
269
270 return kTRUE;
271}
272
273//_____________________________________________________________________
274TH3F* AliTPCCalibKr::CreateHisto(Int_t chamber)
275{
276 //
277 // create new histogram
278 //
279 char name[256];
280 TH3F *h;
281
282 snprintf(name,256,"ADCcluster_ch%d",chamber);
283
284 if( IsIROC(chamber) == kTRUE )
285 {
286 h = new TH3F(name,name,63,0,63,108,0,108,fIrocHistogramNbins,fIrocHistogramMin,fIrocHistogramMax);
287 } else {
288 h = new TH3F(name,name,96,0,96,140,0,140,fOrocHistogramNbins,fOrocHistogramMin,fOrocHistogramMax);
289 }
290 h->SetXTitle("padrow");
291 h->SetYTitle("pad");
292 h->SetZTitle("fADC");
293
294return h;
295}
296
297//_____________________________________________________________________
298Bool_t AliTPCCalibKr::IsIROC(Int_t chamber)
299{
300// check if IROCs
301// returns kTRUE if IROCs and kFALSE if OROCs
302
303 if(chamber>=0 && chamber<36) return kTRUE;
304
305return kFALSE;
306}
307
308//_____________________________________________________________________
309Bool_t AliTPCCalibKr::IsCSide(Int_t chamber)
310{
311// check if C side
312// returns kTRUE if C side and kFALSE if A side
313
314 if((chamber>=18 && chamber<36) || (chamber>=54 && chamber<72)) return kTRUE;
315
316return kFALSE;
317}
318
319//_____________________________________________________________________
320Bool_t AliTPCCalibKr::Update(AliTPCclusterKr *cl)
321{
322 //
323 // fill existing histograms
324 //
325
326 if (!Accept(cl)) return kFALSE;
327 TH3F *h = (TH3F*)fHistoKrArray.At(cl->GetSec());
328 if(!h) return kFALSE;
329
330 h->Fill(cl->GetMax().GetRow(),cl->GetMax().GetPad(),cl->GetADCcluster());
331
332 return kTRUE;
333}
334
335//_____________________________________________________________________
336Bool_t AliTPCCalibKr::Accept(AliTPCclusterKr *cl){
337 //
338 // cuts
339 //
340 /*
341 TCut cutR0("cutR0","fADCcluster/fSize<200"); // adjust it according v seetings -
342 TCut cutR1("cutR1","fADCcluster/fSize>7"); // cosmic tracks and noise removal
343 TCut cutR2("cutR2","fMax.fAdc/fADCcluster<0.4"); // digital noise removal
344 TCut cutR3("cutR3","fMax.fAdc/fADCcluster>0.01"); // noise removal
345 TCut cutR4("cutR4","fMax.fTime>200"); // noise removal
346 TCut cutR5("cutR5","fMax.fTime<600"); // noise removal
347 TCut cutS1("cutS1","fSize<200"); // adjust it according v seetings - cosmic tracks
348 TCut cutAll = cutR0+cutR1+cutR2+cutR3+cutR4+cutR5+cutS1;
349 */
350 /*
351 //R0
352 if ((float)cl->GetADCcluster()/ cl->GetSize() >200) return kFALSE;
353 // R1
354 if ((float)cl->GetADCcluster()/ cl->GetSize() <7) return kFALSE;
355 //R2
356 if ((float)cl->GetMax().GetAdc()/ cl->GetADCcluster() >0.4) return kFALSE;
357 //R3
358 if ((float)cl->GetMax().GetAdc()/ cl->GetADCcluster() <0.01) return kFALSE;
359 //R4
360 if (cl->GetMax().GetTime() < 200) return kFALSE;
361 //R5
362 if (cl->GetMax().GetTime() > 600) return kFALSE;
363 //S1
364 if (cl->GetSize()>200) return kFALSE;
365 if (cl->GetSize()<6) return kFALSE;
366
367 SetADCOverClustSizeRange(7,200);
368 SetMaxADCOverClustADCRange(0.01,0.4);
369 SetTimeRange(200,600);
370 SetClustSizeRange(6,200);
371 */
372
373 //R0
374 if ((float)cl->GetADCcluster()/ cl->GetSize() >fADCOverClustSizeMax) return kFALSE;
375 // R1
376 if ((float)cl->GetADCcluster()/ cl->GetSize() <fADCOverClustSizeMin) return kFALSE;
377 //R2
378 if ((float)cl->GetMax().GetAdc()/ cl->GetADCcluster() >fMaxADCOverClustADCMax) return kFALSE;
379 //R3
380 if ((float)cl->GetMax().GetAdc()/ cl->GetADCcluster() <fMaxADCOverClustADCMin) return kFALSE;
381 //R4
382 if (cl->GetMax().GetTime() > fTimeMax) return kFALSE;
383 //R5
384 if (cl->GetMax().GetTime() < fTimeMin) return kFALSE;
385 //S1
386 if (cl->GetSize()>fClustSizeMax) return kFALSE;
387 if (cl->GetSize()<fClustSizeMin) return kFALSE;
388
389 //
390 // cuts by Adam
391 //
392 /*
393 TCut cutAI0("cutAI0","fTimebinRMS>1.6");
394 TCut cutAI1("cutAI1","fPadRMS>0.825");
395 TCut cutAI2("cutAI2","fRowRMS>0.1");
396 TCut cutAI3("cutAI3","fPads1D<15");
397 TCut cutAI4("cutAI4","fTimebinRMS+11.9-2.15*TMath::Log(1.*fADCcluster)<0");
398
399 TCut cutAIAll = cutAI0+cutAI1+cutAI2+cutAI3+cutAI4;
400
401 TCut cutAO0("cutAO0","fTimebinRMS>0.8");
402 TCut cutAO1("cutAO1","fPadRMS>0.55");
403 TCut cutAO2("cutAO2","fRowRMS>0.1");
404 TCut cutAO3("cutAO3","fPads1D<11");
405 TCut cutAO4("cutAO4","fTimebinRMS+11.9-2.15*TMath::Log(2.*fADCcluster)<0");
406
407 TCut cutAOAll = cutAO0+cutAO1+cutAO2+cutAO3+cutAO4;
408 TCut cutAll("cutAll","(fSec<36&&cutAIAll)||(fSec>=36&&cutAOAll)");
409
410 */
411 if(cl->GetSec()<36){ //IROCs
412 //AI0
413 if((float)cl->GetTimebinRMS() <= fTimebinRmsIrocMin) return kFALSE;
414 //AI1
415 if((float)cl->GetPadRMS() <= fPadRmsIrocMin) return kFALSE;
416 //AI2
417 if((float)cl->GetRowRMS() <= fRowRmsIrocMin) return kFALSE;
418 //AI3
419 if(cl->GetPads1D() >= fClusterPadSize1DIrocMax) return kFALSE;
420 //AI4
421 if((float)cl->GetTimebinRMS()+11.9-2.15*TMath::Log(fCurveCoefficientIroc*(float)cl->GetADCcluster()) >= 0) return kFALSE;
422 }else{ //OROCs
423 //AO0
424 if((float)cl->GetTimebinRMS() <= fTimebinRmsOrocMin) return kFALSE;
425 //AO1
426 if((float)cl->GetPadRMS() <= fPadRmsOrocMin) return kFALSE;
427 //AO2
428 if((float)cl->GetRowRMS() <= fRowRmsOrocMin) return kFALSE;
429 //AO3
430 if(cl->GetPads1D() >= fClusterPadSize1DOrocMax) return kFALSE;
431 //AO4
432 if((float)cl->GetTimebinRMS()+11.9-2.15*TMath::Log(fCurveCoefficientOroc*(float)cl->GetADCcluster()) >= 0) return kFALSE;
433 }
434
435 return kTRUE;
436
437}
438
439//_____________________________________________________________________
440TH3F* AliTPCCalibKr::GetHistoKr(Int_t chamber) const
441{
442 // get histograms from fHistoKrArray
443 return (TH3F*) fHistoKrArray.At(chamber);
444}
445
446//_____________________________________________________________________
447void AliTPCCalibKr::Analyse()
448{
449 //
450 // analyse the histograms and extract krypton calibration parameters
451 //
452
453 // AliTPCCalPads that will contain the calibration parameters
454 AliTPCCalPad* spectrMeanCalPad = new AliTPCCalPad("spectrMean", "spectrMean");
455 AliTPCCalPad* spectrRMSCalPad = new AliTPCCalPad("spectrRMS", "spectrRMS");
456 AliTPCCalPad* fitMeanCalPad = new AliTPCCalPad("fitMean", "fitMean");
457 AliTPCCalPad* fitRMSCalPad = new AliTPCCalPad("fitRMS", "fitRMS");
458 AliTPCCalPad* fitNormChi2CalPad = new AliTPCCalPad("fitNormChi2", "fitNormChi2");
459 AliTPCCalPad* entriesCalPad = new AliTPCCalPad("entries", "entries");
460
461 // file stream for debugging purposes
462 TTreeSRedirector* debugStream = new TTreeSRedirector("debugCalibKr.root");
463
464 // if entries in spectrum less than minEntries, then the fit won't be performed
465 Int_t minEntries = 1; //300;
466
467 Double_t windowFrac = 0.12;
468 // the 3d histogram will be projected on the pads given by the following window size
469 // set the numbers to 0 if you want to do a pad-by-pad calibration
470 UInt_t rowRadius = fRowRadius;//4;
471 UInt_t padRadius = fPadRadius;//4;
472
473 // the step size by which pad and row are incremented is given by the following numbers
474 // set them to 1 if you want the finest granularity
475 UInt_t rowStep = fRowStep;//1;//2 // formerly: 2*rowRadius
476 UInt_t padStep = fPadStep;//1;//4 // formerly: 2*padRadius
477
478 for (Int_t chamber = 0; chamber < 72; chamber++) {
479 //if (chamber != 71) continue;
480 AliTPCCalROC roc(chamber); // I need this only for GetNrows() and GetNPads()
481
482 // Usually I would traverse each pad, take the spectrum for its neighbourhood and
483 // obtain the calibration parameters. This takes very long, so instead I assign the same
484 // calibration values to the whole neighbourhood and then go on to the next neighbourhood.
485 UInt_t nRows = roc.GetNrows();
486 for (UInt_t iRow = 0; iRow < nRows; iRow += rowStep) {
487 UInt_t nPads = roc.GetNPads(iRow);
488 //if (iRow >= 10) break;
489 for (UInt_t iPad = 0; iPad < nPads; iPad += padStep) {
490 //if (iPad >= 20) break;
491 TH3F* h = GetHistoKr(chamber);
492 if (!h) continue;
493
494 // the 3d histogram will be projected on the pads given by the following bounds
495 // for rows and pads
496 Int_t rowLow = iRow - rowRadius;
497 UInt_t rowUp = iRow + rowRadius + rowStep-1;
498 Int_t padLow = iPad - padRadius;
499 UInt_t padUp = iPad + padRadius + padStep-1;
500 // if window goes out of chamber
501 if (rowLow < 0) rowLow = 0;
502 if (rowUp >= nRows) rowUp = nRows - 1;
503 if (padLow < 0) padLow = 0;
504 if (padUp >= nPads) padUp = nPads - 1;
505
506 // project the histogram
507 //TH1D* projH = h->ProjectionZ("projH", rowLow+1, rowUp+1, padLow+1, padUp+1); // SLOW
508 TH1D* projH = ProjectHisto(h, "projH", rowLow+1, rowUp+1, padLow+1, padUp+1);
509
510 // get the number of entries in the spectrum
511 Double_t entries = projH->GetEntries();
512 if (entries < minEntries) { delete projH; continue; }
513
514 // get the two calibration parameters mean of spectrum and RMS of spectrum
515 Double_t histMean = projH->GetMean();
516 Double_t histRMS = (histMean != 0) ? projH->GetRMS() / histMean : 0.;
517
518 // find maximum in spectrum to define a range (given by windowFrac) for which a Gauss is fitted
519 Double_t maxEntries = projH->GetBinCenter(projH->GetMaximumBin());
520 Int_t minBin = projH->FindBin((1.-windowFrac) * maxEntries);
521 Int_t maxBin = projH->FindBin((1.+windowFrac) * maxEntries);
522 Double_t integCharge = projH->Integral(minBin,maxBin);
523
524 Int_t fitResult = projH->Fit("gaus", "Q0", "", (1.-windowFrac) * maxEntries, (1.+windowFrac) * maxEntries);
525
526 if (fitResult != 0) {
527 Error("Analyse", "Error while fitting spectrum for chamber %i, rows %i - %i, pads %i - %i, integrated charge %f.", chamber, rowLow, rowUp, padLow, padUp, integCharge);
528 //printf("[ch%i r%i, p%i] entries = %f, maxEntries = %f, %f\n", chamber, iRow, iPad, entries, maxEntries);
529
530 delete projH;
531 continue;
532 }
533
534 // get the two calibration parameters mean of gauss fit and sigma of gauss fit
535 TF1* gausFit = projH->GetFunction("gaus");
536 Double_t fitMean = gausFit->GetParameter(1);
537 Double_t fitRMS = gausFit->GetParameter(2);
538 Int_t numberFitPoints = gausFit->GetNumberFitPoints();
539 if (numberFitPoints == 0) continue;
540 Double_t fitNormChi2 = gausFit->GetChisquare() / numberFitPoints;
541 delete gausFit;
542 delete projH;
543 if (fitMean <= 0) continue;
544 //printf("[ch%i r%i, p%i] entries = %f, maxEntries = %f, fitMean = %f, fitRMS = %f\n", chamber, iRow, iPad, entries, maxEntries, fitMean, fitRMS);
545
546 // write the calibration parameters for each pad that the 3d histogram was projected onto
547 // (with considering the step size) to the CalPads
548 // rowStep (padStep) odd: round down s/2 and fill this # of rows (pads) in both directions
549 // rowStep (padStep) even: fill s/2 rows (pads) in ascending direction, s/2-1 in descending direction
550 for (Int_t r = iRow - (rowStep/2 - (rowStep+1)%2); r <= (Int_t)(iRow + rowStep/2); r++) {
551 if (r < 0 || r >= (Int_t)nRows) continue;
552 UInt_t nPadsR = roc.GetNPads(r);
553 for (Int_t p = iPad - (padStep/2 - (padStep+1)%2); p <= (Int_t)(iPad + padStep/2); p++) {
554 if (p < 0 || p >= (Int_t)nPadsR) continue;
555 spectrMeanCalPad->GetCalROC(chamber)->SetValue(r, p, histMean);
556 spectrRMSCalPad->GetCalROC(chamber)->SetValue(r, p, histRMS);
557 fitMeanCalPad->GetCalROC(chamber)->SetValue(r, p, fitMean);
558 fitRMSCalPad->GetCalROC(chamber)->SetValue(r, p, fitRMS);
559 fitNormChi2CalPad->GetCalROC(chamber)->SetValue(r, p, fitNormChi2);
560 entriesCalPad->GetCalROC(chamber)->SetValue(r, p, entries);
561
562 (*debugStream) << "calibKr" <<
563 "sector=" << chamber << // chamber number
564 "row=" << r << // row number
565 "pad=" << p << // pad number
566 "histMean=" << histMean << // mean of the spectrum
567 "histRMS=" << histRMS << // RMS of the spectrum divided by the mean
568 "fitMean=" << fitMean << // Gauss fitted mean of the 41.6 keV Kr peak
569 "fitRMS=" << fitRMS << // Gauss fitted sigma of the 41.6 keV Kr peak
570 "fitNormChi2" << fitNormChi2 << // normalized chi square of the Gauss fit
571 "entries=" << entries << // number of entries for the spectrum
572 "\n";
573 }
574 }
575 }
576 }
577 }
578
579 TFile f("calibKr.root", "recreate");
580 spectrMeanCalPad->Write();
581 spectrRMSCalPad->Write();
582 fitMeanCalPad->Write();
583 fitRMSCalPad->Write();
584 fitNormChi2CalPad->Write();
585 entriesCalPad->Write();
586 f.Close();
587 delete spectrMeanCalPad;
588 delete spectrRMSCalPad;
589 delete fitMeanCalPad;
590 delete fitRMSCalPad;
591 delete fitNormChi2CalPad;
592 delete entriesCalPad;
593 delete debugStream;
594}
595
596//_____________________________________________________________________
597TH1D* AliTPCCalibKr::ProjectHisto(TH3F* histo3D, const char* name, Int_t xMin, Int_t xMax, Int_t yMin, Int_t yMax)
598{
599 // project the z-axis of a 3d histo to a specific range of the x- and y-axes,
600 // replaces TH3F::ProjectZ() to gain more speed
601
602 TAxis* xAxis = histo3D->GetXaxis();
603 TAxis* yAxis = histo3D->GetYaxis();
604 TAxis* zAxis = histo3D->GetZaxis();
605 Double_t zMinVal = zAxis->GetXmin();
606 Double_t zMaxVal = zAxis->GetXmax();
607
608 Int_t nBinsZ = zAxis->GetNbins();
609 TH1D* projH = new TH1D(name, name, nBinsZ, zMinVal, zMaxVal);
610
611 Int_t nx = xAxis->GetNbins()+2;
612 Int_t ny = yAxis->GetNbins()+2;
613 Int_t bin = 0;
614 Double_t entries = 0.;
615 for (Int_t x = xMin; x <= xMax; x++) {
616 for (Int_t y = yMin; y <= yMax; y++) {
617 for (Int_t z = 0; z <= nBinsZ+1; z++) {
618 bin = x + nx * (y + ny * z);
619 Double_t val = histo3D->GetBinContent(bin);
620 projH->Fill(zAxis->GetBinCenter(z), val);
621 entries += val;
622 }
623 }
624 }
625 projH->SetEntries((Long64_t)entries);
626 return projH;
627}
628
629//_____________________________________________________________________
630Long64_t AliTPCCalibKr::Merge(TCollection* list) {
631// merge function
632//
633
634if (!list)
635return 0;
636
637if (list->IsEmpty())
638return 1;
639
640TIterator* iter = list->MakeIterator();
641TObject* obj = 0;
642
643 iter->Reset();
644 Int_t count=0;
645 while((obj = iter->Next()) != 0)
646 {
647 AliTPCCalibKr* entry = dynamic_cast<AliTPCCalibKr*>(obj);
648 if (entry == 0) continue;
649
650 for(int i=0; i<72; ++i) {
651 if(IsCSide(i) == kTRUE && fCSide == kTRUE) {
652 ((TH3F*)fHistoKrArray.At(i))->Add( ((TH3F*)entry->fHistoKrArray.At(i)) );
653 }
654
655 if(IsCSide(i) == kFALSE && fASide == kTRUE) {
656 ((TH3F*)fHistoKrArray.At(i))->Add( ((TH3F*)entry->fHistoKrArray.At(i)) );
657 }
658 }
659
660 count++;
661 }
662
663return count;
664}