Line data Source code
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 : /// This class constructs Digits out of Hits
21 : ///
22 : ///
23 :
24 : // --- Standard library ---
25 :
26 : // --- ROOT system ---
27 : #include <TMath.h>
28 : #include <TTree.h>
29 : #include <TMap.h>
30 : #include <TGeoManager.h>
31 : #include <TGeoPhysicalNode.h>
32 : #include <AliGeomManager.h>
33 : #include <TRandom.h>
34 : #include <TF1.h>
35 : #include <TH1F.h>
36 :
37 : // --- AliRoot header files ---
38 : #include "AliRun.h"
39 : #include "AliVZERO.h"
40 : #include "AliVZEROhit.h"
41 : #include "AliRunLoader.h"
42 : #include "AliLoader.h"
43 : #include "AliGRPObject.h"
44 : #include "AliDigitizationInput.h"
45 : #include "AliCDBManager.h"
46 : #include "AliCDBStorage.h"
47 : #include "AliCDBEntry.h"
48 : #include "AliVZEROCalibData.h"
49 : #include "AliCTPTimeParams.h"
50 : #include "AliLHCClockPhase.h"
51 : #include "AliVZEROdigit.h"
52 : #include "AliVZERODigitizer.h"
53 : #include "AliVZEROSDigit.h"
54 : #include "AliLog.h"
55 :
56 12 : ClassImp(AliVZERODigitizer)
57 :
58 : AliVZERODigitizer::AliVZERODigitizer()
59 0 : :AliDigitizer(),
60 0 : fCalibData(GetCalibData()),
61 0 : fPhotoCathodeEfficiency(0.18),
62 0 : fNdigits(0),
63 0 : fDigits(0),
64 0 : fSignalShape(NULL),
65 0 : fPMResponse(NULL),
66 0 : fSinglePhESpectrum(NULL),
67 0 : fEvenOrOdd(kFALSE),
68 0 : fTask(kHits2Digits),
69 0 : fVZERO(NULL)
70 0 : {
71 : // default constructor
72 : // Initialize OCDB and containers used in the digitization
73 :
74 0 : Init();
75 0 : }
76 :
77 : //____________________________________________________________________________
78 : AliVZERODigitizer::AliVZERODigitizer(AliVZERO *vzero, DigiTask_t task)
79 1 : :AliDigitizer(),
80 2 : fCalibData(GetCalibData()),
81 1 : fPhotoCathodeEfficiency(0.18),
82 1 : fNdigits(0),
83 1 : fDigits(0),
84 1 : fSignalShape(NULL),
85 1 : fPMResponse(NULL),
86 1 : fSinglePhESpectrum(NULL),
87 1 : fEvenOrOdd(kFALSE),
88 1 : fTask(task),
89 1 : fVZERO(vzero)
90 5 : {
91 : // constructor
92 : // Initialize OCDB and containers used in the digitization
93 :
94 1 : Init();
95 2 : }
96 :
97 : //____________________________________________________________________________
98 : AliVZERODigitizer::AliVZERODigitizer(AliDigitizationInput* digInput)
99 1 : :AliDigitizer(digInput),
100 2 : fCalibData(GetCalibData()),
101 1 : fPhotoCathodeEfficiency(0.18),
102 1 : fNdigits(0),
103 1 : fDigits(0),
104 1 : fSignalShape(NULL),
105 1 : fPMResponse(NULL),
106 1 : fSinglePhESpectrum(NULL),
107 1 : fEvenOrOdd(kFALSE),
108 1 : fTask(kHits2Digits),
109 1 : fVZERO(NULL)
110 5 : {
111 : // constructor
112 : // Initialize OCDB and containers used in the digitization
113 :
114 1 : Init();
115 2 : }
116 :
117 : //____________________________________________________________________________
118 : AliVZERODigitizer::~AliVZERODigitizer()
119 12 : {
120 : // destructor
121 :
122 2 : if (fDigits) {
123 2 : fDigits->Delete();
124 4 : delete fDigits;
125 2 : fDigits=0;
126 2 : }
127 :
128 2 : if (fSignalShape) {
129 4 : delete fSignalShape;
130 2 : fSignalShape = NULL;
131 2 : }
132 2 : if (fPMResponse) {
133 4 : delete fPMResponse;
134 2 : fPMResponse = NULL;
135 2 : }
136 2 : if (fSinglePhESpectrum) {
137 4 : delete fSinglePhESpectrum;
138 2 : fSinglePhESpectrum = NULL;
139 2 : }
140 :
141 260 : for(Int_t i = 0 ; i < 64; ++i) {
142 384 : if (fTime[i]) delete [] fTime[i];
143 : }
144 6 : }
145 :
146 : //_____________________________________________________________________________
147 : Bool_t AliVZERODigitizer::Init()
148 : {
149 : // Initialises the digitizer
150 : // Initialize OCDB and containers used in the digitization
151 :
152 : // check if the digitizer was already initialized
153 7 : if (fSignalShape) return kTRUE;
154 :
155 4 : fSignalShape = new TF1("VZEROSignalShape",this,&AliVZERODigitizer::SignalShape,0,200,6,"AliVZERODigitizer","SignalShape");
156 : // fSignalShape->SetParameters(0,1.57345e1,-4.25603e-1,2.9,6.40982,3.69339e-01);
157 : // fSignalShape->SetParameters(1.34330e+00,1.13007e+02,-4.95705e-01,
158 : // 3.68911e+00,1.01040e+00, 3.94675e-01);
159 2 : fSignalShape->SetParameters(-1.07335e+00,2.16002e+01,-1.26133e-01,
160 : 1.41619e+00,5.50334e-01,3.86111e-01);
161 4 : fPMResponse = new TF1("VZEROPMResponse",this,&AliVZERODigitizer::PMResponse,-kPMRespTime,2.*kPMRespTime,0,"AliVZERODigitizer","PMResponse");
162 4 : fSinglePhESpectrum = new TF1("VZEROSinglePhESpectrum",this,&AliVZERODigitizer::SinglePhESpectrum,0,20,0,"AliVZERODigitizer","SinglePhESpectrum");
163 :
164 : // Now get the CTP L0->L1 delay
165 4 : AliCDBEntry *entry = AliCDBManager::Instance()->Get("GRP/CTP/CTPtiming");
166 2 : if (!entry) AliFatal("CTP timing parameters are not found in OCDB !");
167 2 : AliCTPTimeParams *ctpParams = (AliCTPTimeParams*)entry->GetObject();
168 2 : Float_t l1Delay = (Float_t)ctpParams->GetDelayL1L0()*25.0;
169 :
170 4 : AliCDBEntry *entry1 = AliCDBManager::Instance()->Get("GRP/CTP/TimeAlign");
171 2 : if (!entry1) AliFatal("CTP time-alignment is not found in OCDB !");
172 2 : AliCTPTimeParams *ctpTimeAlign = (AliCTPTimeParams*)entry1->GetObject();
173 2 : l1Delay += ((Float_t)ctpTimeAlign->GetDelayL1L0()*25.0);
174 :
175 4 : AliCDBEntry *entry2 = AliCDBManager::Instance()->Get("VZERO/Calib/TimeDelays");
176 2 : if (!entry2) AliFatal("VZERO time delays are not found in OCDB !");
177 2 : TH1F *delays = (TH1F*)entry2->GetObject();
178 :
179 4 : AliCDBEntry *entry3 = AliCDBManager::Instance()->Get("GRP/Calib/LHCClockPhase");
180 2 : if (!entry3) AliFatal("LHC clock-phase shift is not found in OCDB !");
181 2 : AliLHCClockPhase *phase = (AliLHCClockPhase*)entry3->GetObject();
182 :
183 260 : for(Int_t i = 0 ; i < 64; ++i) {
184 :
185 5632 : for(Int_t j = 0; j < kNClocks; ++j) fAdc[i][j] = 0;
186 128 : fLeadingTime[i] = fTimeWidth[i] = 0;
187 :
188 128 : fPmGain[i] = fCalibData->GetGain(i);
189 :
190 128 : fAdcPedestal[i][0] = fCalibData->GetPedestal(i);
191 128 : fAdcSigma[i][0] = fCalibData->GetSigma(i);
192 128 : fAdcPedestal[i][1] = fCalibData->GetPedestal(i+64);
193 128 : fAdcSigma[i][1] = fCalibData->GetSigma(i+64);
194 :
195 128 : Int_t board = AliVZEROCalibData::GetBoardNumber(i);
196 384 : fNBins[i] = TMath::Nint(((Float_t)(fCalibData->GetMatchWindow(board)+1)*25.0+
197 256 : (Float_t)(kMaxTDCWidth+1)*fCalibData->GetWidthResolution(board))/
198 128 : fCalibData->GetTimeResolution(board));
199 256 : fNBinsLT[i] = TMath::Nint(((Float_t)(fCalibData->GetMatchWindow(board)+1)*25.0)/
200 128 : fCalibData->GetTimeResolution(board));
201 128 : fBinSize[i] = fCalibData->GetTimeResolution(board);
202 384 : fHptdcOffset[i] = (((Float_t)fCalibData->GetRollOver(board)-
203 384 : (Float_t)fCalibData->GetTriggerCountOffset(board))*25.0+
204 256 : fCalibData->GetTimeOffset(i)-
205 256 : l1Delay-
206 256 : phase->GetMeanPhase()+
207 256 : delays->GetBinContent(i+1)+
208 : kV0Offset);
209 384 : fClockOffset[i] = (((Float_t)fCalibData->GetRollOver(board)-
210 384 : (Float_t)fCalibData->GetTriggerCountOffset(board))*25.0+
211 256 : fCalibData->GetTimeOffset(i)-
212 128 : l1Delay+
213 : kV0Offset);
214 : // Needed for Run2 data simulation
215 128 : if (AliCDBManager::Instance()->GetRun() >= 215011) {
216 0 : fHptdcOffset[i] -= 1.4;
217 0 : fClockOffset[i] += 269.;
218 0 : }
219 :
220 128 : fTime[i] = new Float_t[fNBins[i]];
221 128 : memset(fTime[i],0,fNBins[i]*sizeof(Float_t));
222 : }
223 :
224 : return kTRUE;
225 3 : }
226 :
227 : //____________________________________________________________________________
228 : void AliVZERODigitizer::Digitize(Option_t* /*option*/)
229 : {
230 : // Creates digits from hits
231 10 : fNdigits = 0;
232 :
233 6 : if (fVZERO && !fDigInput) {
234 1 : AliLoader *loader = fVZERO->GetLoader();
235 1 : if (!loader) {
236 0 : AliError("Can not get VZERO Loader via AliVZERO object!");
237 0 : return;
238 : }
239 1 : AliRunLoader* runLoader = AliRunLoader::Instance();
240 10 : for (Int_t iEvent = 0; iEvent < runLoader->GetNumberOfEvents(); ++iEvent) {
241 4 : runLoader->GetEvent(iEvent);
242 8 : if (fTask == kHits2Digits) {
243 4 : DigitizeHits();
244 0 : DigitizeSDigits();
245 0 : WriteDigits(loader);
246 0 : }
247 : else {
248 : DigitizeHits();
249 4 : WriteSDigits(loader);
250 : }
251 : }
252 1 : }
253 4 : else if (fDigInput) {
254 4 : ReadSDigits();
255 4 : DigitizeSDigits();
256 4 : AliRunLoader *currentLoader = AliRunLoader::GetRunLoader(fDigInput->GetOutputFolderName());
257 4 : AliLoader *loader = currentLoader->GetLoader("VZEROLoader");
258 4 : if (!loader) {
259 0 : AliError("Cannot get VZERO Loader via RunDigitizer!");
260 0 : return;
261 : }
262 4 : WriteDigits(loader);
263 4 : }
264 : else {
265 0 : AliFatal("Invalid digitization task! Exiting!");
266 : }
267 4 : }
268 :
269 : //____________________________________________________________________________
270 : void AliVZERODigitizer::AddDigit(Int_t pmnumber, Float_t time, Float_t width, Bool_t integrator, Short_t *chargeADC, Int_t *labels)
271 : {
272 :
273 : // Adds Digit
274 :
275 256 : TClonesArray &ldigits = *fDigits;
276 :
277 256 : new(ldigits[fNdigits++]) AliVZEROdigit(pmnumber,time,width,integrator,chargeADC,labels);
278 :
279 256 : }
280 : //____________________________________________________________________________
281 : void AliVZERODigitizer::AddSDigit(Int_t pmnumber, Int_t nbins, Float_t *charges, Int_t *labels)
282 : {
283 :
284 : // Adds SDigit
285 :
286 512 : TClonesArray &ldigits = *fDigits;
287 :
288 256 : new(ldigits[fNdigits++]) AliVZEROSDigit(pmnumber,nbins,charges,labels);
289 :
290 256 : }
291 : //____________________________________________________________________________
292 : void AliVZERODigitizer::ResetDigits()
293 : {
294 :
295 : // Clears Digits
296 :
297 16 : fNdigits = 0;
298 16 : if (fDigits) fDigits->Clear();
299 8 : }
300 :
301 : //____________________________________________________________________________
302 : AliVZEROCalibData* AliVZERODigitizer::GetCalibData() const
303 :
304 : {
305 4 : AliCDBManager *man = AliCDBManager::Instance();
306 :
307 : AliCDBEntry *entry=0;
308 :
309 4 : entry = man->Get("VZERO/Calib/Data");
310 :
311 : // if(!entry){
312 : // AliWarning("Load of calibration data from default storage failed!");
313 : // AliWarning("Calibration data will be loaded from local storage ($ALICE_ROOT)");
314 : // Int_t runNumber = man->GetRun();
315 : // entry = man->GetStorage("local://$ALICE_ROOT/OCDB")
316 : // ->Get("VZERO/Calib/Data",runNumber);
317 : //
318 : // }
319 :
320 : // Retrieval of data in directory VZERO/Calib/Data:
321 :
322 :
323 : AliVZEROCalibData *calibdata = 0;
324 :
325 4 : if (entry) calibdata = (AliVZEROCalibData*) entry->GetObject();
326 2 : if (!calibdata) AliFatal("No calibration data from calibration database !");
327 :
328 2 : return calibdata;
329 :
330 0 : }
331 :
332 : double AliVZERODigitizer::SignalShape(double *x, double *par)
333 : {
334 : // this function simulates the time
335 : // of arrival of the photons at the
336 : // photocathode
337 19840 : Double_t xx = x[0];
338 9920 : if (xx <= par[0]) return 0;
339 9920 : Double_t a = 1./TMath::Power((xx-par[0])/par[1],1./par[2]);
340 10091 : if (xx <= par[3]) return a;
341 9749 : Double_t b = 1./TMath::Power((xx-par[3])/par[4],1./par[5]);
342 9749 : Double_t f = a*b/(a+b);
343 29247 : AliDebug(100,Form("x=%f func=%f",xx,f));
344 : return f;
345 9920 : }
346 :
347 : double AliVZERODigitizer::PMResponse(double *x, double * /* par */)
348 : {
349 : // this function describes the
350 : // PM time response to a single
351 : // photoelectron
352 493322 : Double_t xx = x[0]+kPMRespTime;
353 246661 : return xx*xx*TMath::Exp(-xx*xx/(kPMRespTime*kPMRespTime));
354 : }
355 :
356 : double AliVZERODigitizer::SinglePhESpectrum(double *x, double * /* par */)
357 : {
358 : // this function describes the
359 : // PM amplitude response to a single
360 : // photoelectron
361 11520 : Double_t xx = x[0];
362 5760 : if (xx < 0) return 0;
363 5760 : return (TMath::Poisson(xx,kPMNbOfSecElec)+kPMTransparency*TMath::Poisson(xx,1.0));
364 5760 : }
365 :
366 : Int_t AliVZERODigitizer::Cell2Pmt(Int_t cell) const
367 : {
368 : // The method maps the scintillator
369 : // indexes to the PM ones
370 128 : if (cell < 0 || cell >= 80) {
371 0 : AliError(Form("Wrong VZERO cell index %d",cell));
372 0 : return -1;
373 : }
374 71 : if (cell < 16) return cell;
375 70 : if (cell < 48) return 8 + cell/2;
376 44 : return cell - 16;
377 64 : }
378 :
379 : void AliVZERODigitizer::DigitizeHits()
380 : {
381 : // Digitize the hits to the level of
382 : // SDigits (fTime arrays)
383 :
384 524 : for(Int_t i = 0 ; i < 64; ++i) {
385 256 : memset(fTime[i],0,fNBins[i]*sizeof(Float_t));
386 256 : fLabels[i][0] = fLabels[i][1] = fLabels[i][2] = -1;
387 : }
388 4 : Float_t integral = fPMResponse->Integral(-kPMRespTime,2.*kPMRespTime);
389 4 : Float_t meansPhE = fSinglePhESpectrum->Mean(0,20);
390 :
391 4 : AliLoader* loader = fVZERO->GetLoader();
392 4 : if (!loader) {
393 0 : AliError("Can not get VZERO Loader!");
394 0 : return;
395 : }
396 4 : loader->LoadHits();
397 4 : TTree* treeH = loader->TreeH();
398 4 : if (!treeH) {
399 0 : AliError("Cannot get TreeH!");
400 0 : return;
401 : }
402 4 : TClonesArray* hits = fVZERO->Hits();
403 :
404 : // Now makes Digits from hits
405 4 : Int_t nTracks = (Int_t) treeH->GetEntries();
406 232 : for (Int_t iTrack = 0; iTrack < nTracks; iTrack++) {
407 112 : fVZERO->ResetHits();
408 112 : treeH->GetEvent(iTrack);
409 112 : Int_t nHits = hits->GetEntriesFast();
410 352 : for (Int_t iHit = 0; iHit < nHits; iHit++) {
411 64 : AliVZEROhit* hit = (AliVZEROhit *)hits->UncheckedAt(iHit);
412 64 : Int_t nPhot = hit->Nphot();
413 64 : Int_t cell = hit->Cell();
414 64 : Int_t pmt = Cell2Pmt(cell);
415 64 : if (pmt < 0) continue;
416 64 : Int_t trackLabel = hit->GetTrack();
417 178 : for(Int_t l = 0; l < 3; ++l) {
418 86 : if (fLabels[pmt][l] < 0) {
419 62 : fLabels[pmt][l] = trackLabel;
420 62 : break;
421 : }
422 : }
423 : Float_t dt_scintillator;
424 128 : if (AliCDBManager::Instance()->GetRun() < 215011) {
425 128 : dt_scintillator = gRandom->Gaus(0,kIntTimeRes);
426 64 : }
427 : else {
428 0 : dt_scintillator = gRandom->Gaus(0,kIntTimeResRing[pmt/8]);
429 : }
430 64 : Float_t t = dt_scintillator + 1e9*hit->Tof();
431 84 : if (pmt < 32) t += kV0CDelayCables;
432 64 : t += fHptdcOffset[pmt];
433 : Int_t nPhE;
434 64 : Float_t prob = fCalibData->GetLightYields(pmt)*fPhotoCathodeEfficiency; // Optical losses included!
435 128 : if (nPhot > 100)
436 192 : nPhE = (Int_t)gRandom->Gaus(prob*Float_t(nPhot)+0.5,
437 64 : sqrt(Float_t(nPhot)*prob*(1.-prob)));
438 : else
439 0 : nPhE = gRandom->Binomial(nPhot,prob);
440 64 : Float_t charge = TMath::Qe()*fPmGain[pmt]*fBinSize[pmt]/integral;
441 2804 : for (Int_t iPhE = 0; iPhE < nPhE; ++iPhE) {
442 1338 : Float_t tPhE = t + fSignalShape->GetRandom(0,fBinSize[pmt]*Float_t(fNBins[pmt]));
443 1338 : Float_t gainVar = fSinglePhESpectrum->GetRandom(0,20)/meansPhE;
444 1338 : Int_t firstBin = TMath::Max(0,(Int_t)((tPhE-kPMRespTime)/fBinSize[pmt]));
445 1338 : Int_t lastBin = TMath::Min(fNBins[pmt]-1,(Int_t)((tPhE+2.*kPMRespTime)/fBinSize[pmt]));
446 494846 : for(Int_t iBin = firstBin; iBin <= lastBin; ++iBin) {
447 246085 : Float_t tempT = fBinSize[pmt]*(0.5+iBin)-tPhE;
448 246085 : fTime[pmt][iBin] += gainVar*charge*fPMResponse->Eval(tempT);
449 : }
450 : } // ph.e. loop
451 64 : } // hit loop
452 : } // track loop
453 4 : loader->UnloadHits();
454 8 : }
455 :
456 :
457 : void AliVZERODigitizer::DigitizeSDigits()
458 : {
459 : // Digitize the fTime arrays (SDigits) to the level of
460 : // Digits (fAdc arrays)
461 524 : for(Int_t i = 0 ; i < 64; ++i) {
462 11264 : for(Int_t j = 0; j < kNClocks; ++j) fAdc[i][j] = 0;
463 256 : fLeadingTime[i] = fTimeWidth[i] = 0;
464 : }
465 :
466 4 : Float_t maximum = 0.9*fSignalShape->GetMaximum(0,200); // Not exact, one needs to do this on the convoluted
467 4 : Float_t integral2 = fSignalShape->Integral(0,200); // function. Anyway the effect is small <10% on the 2.5 ADC thr
468 520 : for (Int_t ipmt = 0; ipmt < 64; ++ipmt) {
469 256 : Float_t thr = fCalibData->GetCalibDiscriThr(ipmt,kFALSE,AliCDBManager::Instance()->GetRun())*kChargePerADC*maximum*fBinSize[ipmt]/integral2;
470 : Bool_t ltFound = kFALSE, ttFound = kFALSE;
471 918016 : for (Int_t iBin = 0; iBin < fNBins[ipmt]; ++iBin) {
472 458752 : Float_t t = fBinSize[ipmt]*Float_t(iBin);
473 917504 : if (fTime[ipmt][iBin] > thr) {
474 465524 : if (!ltFound && (iBin < fNBinsLT[ipmt])) {
475 : ltFound = kTRUE;
476 41 : fLeadingTime[ipmt] = t;
477 41 : }
478 : }
479 : else {
480 452021 : if (ltFound) {
481 31914 : if (!ttFound) {
482 : ttFound = kTRUE;
483 41 : fTimeWidth[ipmt] = t - fLeadingTime[ipmt];
484 41 : }
485 : }
486 : }
487 458752 : Float_t tadc = t - fClockOffset[ipmt];
488 458752 : Int_t clock = kNClocks/2 - Int_t(tadc/25.0);
489 458752 : if (clock >= 0 && clock < kNClocks)
490 458752 : fAdc[ipmt][clock] += fTime[ipmt][iBin]/kChargePerADC;
491 : }
492 768 : AliDebug(1,Form("Channel %d Offset %f Time %f",ipmt,fClockOffset[ipmt],fLeadingTime[ipmt]));
493 256 : Int_t board = AliVZEROCalibData::GetBoardNumber(ipmt);
494 297 : if (ltFound && ttFound) {
495 82 : fTimeWidth[ipmt] = fCalibData->GetWidthResolution(board)*
496 41 : Float_t(Int_t(fTimeWidth[ipmt]/fCalibData->GetWidthResolution(board)));
497 41 : if (fTimeWidth[ipmt] < Float_t(kMinTDCWidth)*fCalibData->GetWidthResolution(board))
498 0 : fTimeWidth[ipmt] = Float_t(kMinTDCWidth)*fCalibData->GetWidthResolution(board);
499 41 : if (fTimeWidth[ipmt] > Float_t(kMaxTDCWidth)*fCalibData->GetWidthResolution(board))
500 0 : fTimeWidth[ipmt] = Float_t(kMaxTDCWidth)*fCalibData->GetWidthResolution(board);
501 : }
502 : }
503 :
504 4 : fEvenOrOdd = gRandom->Integer(2);
505 520 : for (Int_t j=0; j<64; ++j){
506 11264 : for (Int_t iClock = 0; iClock < kNClocks; ++iClock) {
507 5376 : Int_t integrator = (iClock + fEvenOrOdd) % 2;
508 16128 : AliDebug(1,Form("ADC %d %d %f",j,iClock,fAdc[j][iClock]));
509 5376 : fAdc[j][iClock] += gRandom->Gaus(fAdcPedestal[j][integrator], fAdcSigma[j][integrator]);
510 : }
511 : }
512 :
513 4 : }
514 :
515 : void AliVZERODigitizer::WriteDigits(AliLoader *loader)
516 : {
517 : // Take fAdc arrays filled by the previous
518 : // method and produce and add digits to the digit Tree
519 :
520 8 : loader->LoadDigits("UPDATE");
521 :
522 8 : if (!loader->TreeD()) loader->MakeTree("D");
523 4 : loader->MakeDigitsContainer();
524 4 : TTree* treeD = loader->TreeD();
525 4 : DigitsArray();
526 4 : if (treeD->GetBranch("VZERODigits"))
527 0 : treeD->SetBranchAddress("VZERODigits", &fDigits);
528 : else
529 4 : treeD->Bronch("VZERODigit", "TClonesArray", &fDigits);
530 :
531 4 : Short_t *chargeADC = new Short_t[kNClocks];
532 520 : for (Int_t i=0; i<64; i++) {
533 11264 : for (Int_t j = 0; j < kNClocks; ++j) {
534 5376 : Int_t tempadc = Int_t(fAdc[i][j]);
535 5376 : if (tempadc > 1023) tempadc = 1023;
536 5376 : chargeADC[j] = tempadc;
537 : }
538 256 : AddDigit(i, fLeadingTime[i], fTimeWidth[i], Bool_t((10+fEvenOrOdd)%2), chargeADC, fLabels[i]);
539 : }
540 8 : delete [] chargeADC;
541 :
542 4 : treeD->Fill();
543 4 : loader->WriteDigits("OVERWRITE");
544 4 : loader->UnloadDigits();
545 4 : ResetDigits();
546 4 : }
547 :
548 : void AliVZERODigitizer::WriteSDigits(AliLoader *loader)
549 : {
550 : // Take fTime arrays filled by the previous
551 : // method and produce and add sdigits to the sdigit Tree
552 :
553 8 : loader->LoadSDigits("UPDATE");
554 :
555 8 : if (!loader->TreeS()) loader->MakeTree("S");
556 4 : loader->MakeSDigitsContainer();
557 4 : TTree* treeS = loader->TreeS();
558 4 : SDigitsArray();
559 4 : if (treeS->GetBranch("VZEROSDigits"))
560 0 : treeS->SetBranchAddress("VZEROSDigits", &fDigits);
561 : else
562 4 : treeS->Bronch("VZEROSDigit", "TClonesArray", &fDigits);
563 :
564 520 : for (Int_t ipmt = 0; ipmt < 64; ++ipmt) {
565 256 : AddSDigit(ipmt,fNBins[ipmt],fTime[ipmt],fLabels[ipmt]);
566 : }
567 :
568 4 : treeS->Fill();
569 4 : loader->WriteSDigits("OVERWRITE");
570 4 : loader->UnloadSDigits();
571 4 : ResetDigits();
572 4 : }
573 :
574 : void AliVZERODigitizer::ReadSDigits()
575 : {
576 : // Read SDigits which are then to precessed
577 : // in the following method
578 524 : for(Int_t i = 0 ; i < 64; ++i) {
579 256 : memset(fTime[i],0,fNBins[i]*sizeof(Float_t));
580 256 : fLabels[i][0] = fLabels[i][1] = fLabels[i][2] = -1;
581 : }
582 :
583 : // Loop over input files
584 4 : Int_t nFiles= fDigInput->GetNinputs();
585 20 : for (Int_t inputFile = 0; inputFile < nFiles; inputFile++) {
586 : // Get the current loader
587 : AliRunLoader* currentLoader =
588 4 : AliRunLoader::GetRunLoader(fDigInput->GetInputFolderName(inputFile));
589 :
590 4 : AliLoader *loader = currentLoader->GetLoader("VZEROLoader");
591 4 : loader->LoadSDigits("READ");
592 :
593 : // Get the tree of summable digits
594 4 : TTree* sdigitsTree = loader->TreeS();
595 4 : if (!sdigitsTree) {
596 0 : AliError("No sdigit tree from digInput");
597 0 : continue;
598 : }
599 :
600 : // Get the branch
601 4 : TBranch* sdigitsBranch = sdigitsTree->GetBranch("VZEROSDigit");
602 4 : if (!sdigitsBranch) {
603 0 : AliError("Failed to get sdigit branch");
604 0 : return;
605 : }
606 :
607 : // Set the branch address
608 4 : TClonesArray *sdigitsArray = NULL;
609 4 : sdigitsBranch->SetAddress(&sdigitsArray);
610 :
611 : // Sum contributions from the sdigits
612 : // Get number of entries in the tree
613 4 : Int_t nentries = Int_t(sdigitsBranch->GetEntries());
614 16 : for (Int_t entry = 0; entry < nentries; ++entry) {
615 4 : sdigitsBranch->GetEntry(entry);
616 : // Get the number of sdigits
617 4 : Int_t nsdigits = sdigitsArray->GetEntries();
618 520 : for (Int_t sdigit = 0; sdigit < nsdigits; sdigit++) {
619 256 : AliVZEROSDigit* sDigit = static_cast<AliVZEROSDigit*>(sdigitsArray->UncheckedAt(sdigit));
620 256 : Int_t pmNumber = sDigit->PMNumber();
621 256 : Int_t nbins = sDigit->GetNBins();
622 256 : if (nbins != fNBins[pmNumber]) {
623 0 : AliError(Form("Incompatible number of bins between digitizer (%d) and sdigit (%d) for PM %d! Skipping sdigit!",
624 : fNBins[pmNumber],nbins,pmNumber));
625 0 : continue;
626 : }
627 : // Sum the charges
628 256 : Float_t *charges = sDigit->GetCharges();
629 918016 : for(Int_t iBin = 0; iBin < nbins; ++iBin) fTime[pmNumber][iBin] += charges[iBin];
630 : // and the labels
631 256 : Int_t *labels = sDigit->GetTracks();
632 : Int_t j = 0;
633 638 : for(Int_t i = 0; i < 3; ++i) {
634 316 : if (fLabels[pmNumber][i] < 0) {
635 570 : if (labels[j] < 0) break;
636 62 : fLabels[pmNumber][i] = labels[j];
637 62 : j++;
638 62 : }
639 : }
640 256 : }
641 : }
642 4 : loader->UnloadSDigits();
643 4 : }
644 8 : }
645 :
646 : //____________________________________________________________________
647 : TClonesArray*
648 : AliVZERODigitizer::DigitsArray()
649 : {
650 : // Initialize digit array if not already and
651 : // return pointer to it.
652 8 : if (!fDigits) {
653 2 : fDigits = new TClonesArray("AliVZEROdigit", 64);
654 1 : fNdigits = 0;
655 1 : }
656 4 : return fDigits;
657 0 : }
658 :
659 : //____________________________________________________________________
660 : TClonesArray*
661 : AliVZERODigitizer::SDigitsArray()
662 : {
663 : // Initialize sdigit array if not already and
664 : // return pointer to it.
665 8 : if (!fDigits) {
666 2 : fDigits = new TClonesArray("AliVZEROSDigit", 64);
667 1 : fNdigits = 0;
668 1 : }
669 4 : return fDigits;
670 0 : }
|
