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