]> git.uio.no Git - u/mrichter/AliRoot.git/blob - ZDC/AliZDCReconstructor.cxx
ZEM signal not yet splitted in ZDC ESD
[u/mrichter/AliRoot.git] / ZDC / AliZDCReconstructor.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$ */
17
18 ///////////////////////////////////////////////////////////////////////////////
19 //                                                                           //
20 // class for ZDC reconstruction                                              //
21 //                                                                           //
22 ///////////////////////////////////////////////////////////////////////////////
23
24
25 #include <TF1.h>
26
27 #include "AliRunLoader.h"
28 #include "AliRawReader.h"
29 #include "AliESDEvent.h"
30 #include "AliZDCDigit.h"
31 #include "AliZDCRawStream.h"
32 #include "AliZDCReco.h"
33 #include "AliZDCReconstructor.h"
34 #include "AliZDCCalibData.h"
35
36
37 ClassImp(AliZDCReconstructor)
38
39
40 //_____________________________________________________________________________
41 AliZDCReconstructor:: AliZDCReconstructor() :
42
43   fZNCen(new TF1("fZNCen", 
44         "(-2.287920+sqrt(2.287920*2.287920-4*(-0.007629)*(11.921710-x)))/(2*(-0.007629))",0.,164.)),
45   fZNPer(new TF1("fZNPer",
46       "(-37.812280-sqrt(37.812280*37.812280-4*(-0.190932)*(-1709.249672-x)))/(2*(-0.190932))",0.,164.)),
47   fZPCen(new TF1("fZPCen",
48        "(-1.321353+sqrt(1.321353*1.321353-4*(-0.007283)*(3.550697-x)))/(2*(-0.007283))",0.,60.)),
49   fZPPer(new TF1("fZPPer",
50       "(-42.643308-sqrt(42.643308*42.643308-4*(-0.310786)*(-1402.945615-x)))/(2*(-0.310786))",0.,60.)),
51   fZDCCen(new TF1("fZDCCen",
52       "(-1.934991+sqrt(1.934991*1.934991-4*(-0.004080)*(15.111124-x)))/(2*(-0.004080))",0.,225.)),
53   fZDCPer(new TF1("fZDCPer",
54       "(-34.380639-sqrt(34.380639*34.380639-4*(-0.104251)*(-2612.189017-x)))/(2*(-0.104251))",0.,225.)),
55   fbCen(new TF1("fbCen","-0.056923+0.079703*x-0.0004301*x*x+0.000001366*x*x*x",0.,220.)),
56   fbPer(new TF1("fbPer","17.943998-0.046846*x+0.000074*x*x",0.,220.)),
57   //
58   fZEMn(new TF1("fZEMn","121.7-0.1934*x+0.00007565*x*x",0.,1200.)),
59   fZEMp(new TF1("fZEMp","80.05-0.1315*x+0.00005327*x*x",0.,1200.)),
60   fZEMsp(new TF1("fZEMsp","201.7-0.325*x+0.0001292*x*x",0.,1200.)),
61   fZEMb(new TF1("fZEMb",
62         "13.83-0.02851*x+5.101e-5*x*x-7.305e-8*x*x*x+5.101e-11*x*x*x*x-1.25e-14*x*x*x*x*x",0.,1200.)),
63   //
64   fCalibData(GetCalibData())
65
66 {
67   // **** Default constructor
68
69 }
70
71
72 //_____________________________________________________________________________
73 AliZDCReconstructor::~AliZDCReconstructor()
74 {
75 // destructor
76
77   delete fZNCen;
78   delete fZNPer;
79   delete fZPCen;
80   delete fZPPer;
81   delete fZDCCen;
82   delete fZDCPer;
83   delete fbCen;
84   delete fbPer;
85   delete fZEMn;
86   delete fZEMp;
87   delete fZEMsp;
88   delete fZEMb;
89
90 }
91
92
93 //_____________________________________________________________________________
94 void AliZDCReconstructor::Reconstruct(TTree* digitsTree, TTree* clustersTree) const
95 {
96   // *** Local ZDC reconstruction for digits
97   // Works on the current event
98     
99   // Retrieving calibration data  
100   Float_t meanPed[47];
101   for(Int_t jj=0; jj<47; jj++) meanPed[jj] = fCalibData->GetMeanPed(jj);
102
103   // get digits
104   AliZDCDigit digit;
105   AliZDCDigit* pdigit = &digit;
106   digitsTree->SetBranchAddress("ZDC", &pdigit);
107
108   // loop over digits
109   Float_t tZN1CorrHG[]={0.,0.,0.,0.,0.}, tZP1CorrHG[]={0.,0.,0.,0.,0.}; 
110   Float_t dZEMCorrHG=0.; 
111   Float_t tZN2CorrHG[]={0.,0.,0.,0.,0.}, tZP2CorrHG[]={0.,0.,0.,0.,0.};
112   Float_t tZN1CorrLG[]={0.,0.,0.,0.,0.}, tZP1CorrLG[]={0.,0.,0.,0.,0.};
113   Float_t dZEMCorrLG=0.; 
114   Float_t tZN2CorrLG[]={0.,0.,0.,0.,0.}, tZP2CorrLG[]={0.,0.,0.,0.,0.};
115   
116   //printf("\n\t # of digits in tree: %d\n",(Int_t) digitsTree->GetEntries());
117   for (Int_t iDigit = 0; iDigit < (digitsTree->GetEntries()/2); iDigit++) {
118     digitsTree->GetEntry(iDigit);
119     if (!pdigit) continue;
120     //pdigit->Print("");
121     //  
122     Int_t det = digit.GetSector(0);
123     Int_t quad = digit.GetSector(1);
124     Int_t pedindex = -1;
125     //printf("\n\t #%d det %d quad %d", iDigit, det, quad);
126     //
127     if(det == 1){ // *** ZN1
128        pedindex = quad;
129        tZN1CorrHG[quad] = (Float_t) (digit.GetADCValue(0)-meanPed[pedindex]);
130        if(tZN1CorrHG[quad]<0.) tZN1CorrHG[quad] = 0.;
131        tZN1CorrLG[quad] = (Float_t) (digit.GetADCValue(1)-meanPed[pedindex+5]);
132        if(tZN1CorrLG[quad]<0.) tZN1CorrLG[quad] = 0.;
133        //printf("\t pedindex %d tZN1CorrHG[%d] = %1.0f tZN1CorrLG[%d] = %1.0f", 
134        //       pedindex, quad, tZN1CorrHG[quad], quad, tZN1CorrLG[quad]);
135     }
136     else if(det == 2){ // *** ZP1
137        pedindex = quad+10;
138        tZP1CorrHG[quad] = (Float_t) (digit.GetADCValue(0)-meanPed[pedindex]);
139        if(tZP1CorrLG[quad]<0.) tZP1CorrLG[quad] = 0.;
140        tZP1CorrLG[quad] = (Float_t) (digit.GetADCValue(1)-meanPed[pedindex+5]);
141        if(tZP1CorrHG[quad]<0.) tZP1CorrHG[quad] = 0.;
142        //printf("\t pedindex %d tZP1CorrHG[%d] = %1.0f tZP1CorrLG[%d] = %1.0f", 
143        //       pedindex, quad, tZP1CorrHG[quad], quad, tZP1CorrLG[quad]);
144     }
145     else if(det == 3){
146        if(quad == 1){       // *** ZEM1  
147          pedindex = quad+19;
148          dZEMCorrHG += (Float_t) (digit.GetADCValue(0)-meanPed[pedindex]); 
149          if(dZEMCorrHG<0.) dZEMCorrHG = 0.;
150          dZEMCorrLG += (Float_t) (digit.GetADCValue(1)-meanPed[pedindex+2]); 
151          if(dZEMCorrLG<0.) dZEMCorrLG = 0.;
152          //printf("\t pedindex %d ADC(0) = %d ped = %1.0f ADCCorr = %1.0f\n", 
153          //     pedindex, digit.GetADCValue(0), meanPed[pedindex], dZEMCorrHG);
154          //printf("\t pedindex %d ADC(1) = %d ped = %1.0f ADCCorr = %1.0f\n", 
155          //     pedindex+2, digit.GetADCValue(1), meanPed[pedindex+2], dZEMCorrLG);
156          ////printf("\t pedindex %d dZEMCorrHG = %1.0f dZEMCorrLG = %1.0f\n", pedindex, dZEMCorrHG, dZEMCorrLG);
157        }
158        else if(quad == 2){  // *** ZEM1
159          pedindex = quad+19;
160          dZEMCorrHG += (Float_t) (digit.GetADCValue(0)-meanPed[pedindex]); 
161          if(dZEMCorrHG<0.) dZEMCorrHG = 0.;
162          dZEMCorrLG += (Float_t) (digit.GetADCValue(1)-meanPed[pedindex+2]); 
163          if(dZEMCorrLG<0.) dZEMCorrLG = 0.;
164          //printf("\t pedindex %d ADC(0) = %d ped = %1.0f ADCCorr = %1.0f\n", 
165          //     pedindex, digit.GetADCValue(0), meanPed[pedindex], dZEMCorrHG);
166          //printf("\t pedindex %d ADC(1) = %d ped = %1.0f ADCCorr = %1.0f\n", 
167          //     pedindex+2, digit.GetADCValue(1),meanPed[pedindex+2], dZEMCorrLG);
168          ////printf("\t pedindex %d dZEMCorrHG = %1.0f dZEMCorrLG = %1.0f\n", pedindex, dZEMCorrHG, dZEMCorrLG);
169        }
170     }
171     else if(det == 4){  // *** ZN2
172        pedindex = quad+24;
173        tZN2CorrHG[quad] = (Float_t) (digit.GetADCValue(0)-meanPed[pedindex]);
174        if(tZN2CorrHG[quad]<0.) tZN2CorrHG[quad] = 0.;
175        tZN2CorrLG[quad] = (Float_t) (digit.GetADCValue(1)-meanPed[pedindex+5]);
176        if(tZN2CorrLG[quad]<0.) tZN2CorrLG[quad] = 0.;
177        //printf("\t pedindex %d tZN2CorrHG[%d] = %1.0f tZN2CorrLG[%d] = %1.0f\n", 
178        //       pedindex, quad, tZN2CorrHG[quad], quad, tZN2CorrLG[quad]);
179     }
180     else if(det == 5){  // *** ZP2 
181        pedindex = quad+34;
182        tZP2CorrHG[quad] = (Float_t) (digit.GetADCValue(0)-meanPed[pedindex]);
183        if(tZP2CorrHG[quad]<0.) tZP2CorrHG[quad] = 0.;
184        tZP2CorrLG[quad] = (Float_t) (digit.GetADCValue(1)-meanPed[pedindex+5]);
185        if(tZP2CorrLG[quad]<0.) tZP2CorrLG[quad] = 0.;
186        //printf("\t pedindex %d tZP2CorrHG[%d] = %1.0f tZP2CorrLG[%d] = %1.0f\n", 
187        //       pedindex, quad, tZP2CorrHG[quad], quad, tZP2CorrLG[quad]);
188     }
189   }
190
191   // reconstruct the event
192     ReconstructEvent(clustersTree, tZN1CorrHG, tZP1CorrHG, tZN2CorrHG, 
193         tZP2CorrHG, tZN1CorrLG, tZP1CorrLG, tZN2CorrLG, 
194         tZP2CorrLG, dZEMCorrHG);
195
196 }
197
198 //_____________________________________________________________________________
199 void AliZDCReconstructor::Reconstruct(AliRawReader* rawReader, TTree* clustersTree) const
200 {
201   // *** ZDC raw data reconstruction
202   // Works on the current event
203   
204   // Retrieving calibration data  
205   Float_t meanPed[47];
206   for(Int_t jj=0; jj<47; jj++) meanPed[jj] = fCalibData->GetMeanPed(jj);
207
208   rawReader->Reset();
209
210   // loop over raw data rawDatas
211   Float_t tZN1CorrHG[]={0.,0.,0.,0.,0.}, tZP1CorrHG[]={0.,0.,0.,0.,0.};
212   Float_t dZEMCorrHG=0.;
213   Float_t tZN2CorrHG[]={0.,0.,0.,0.,0.}, tZP2CorrHG[]={0.,0.,0.,0.,0.};
214   Float_t tZN1CorrLG[]={0.,0.,0.,0.,0.}, tZP1CorrLG[]={0.,0.,0.,0.,0.};
215   Float_t dZEMCorrLG=0.; 
216   Float_t tZN2CorrLG[]={0.,0.,0.,0.,0.}, tZP2CorrLG[]={0.,0.,0.,0.,0.};
217   //
218   AliZDCRawStream rawData(rawReader);
219   while (rawData.Next()) {
220     if(rawData.IsADCDataWord()){
221       Int_t det = rawData.GetSector(0);
222       Int_t quad = rawData.GetSector(1);
223       Int_t gain = rawData.GetADCGain();
224       Int_t pedindex;
225       //
226       if(det == 1){    
227         pedindex = quad;
228         if(gain == 0) tZN1CorrHG[quad]  += (Float_t) (rawData.GetADCValue()-meanPed[pedindex]); 
229         else tZN1CorrLG[quad]  += (Float_t) (rawData.GetADCValue()-meanPed[pedindex+5]); 
230       }
231       else if(det == 2){ 
232         pedindex = quad+10;
233         if(gain == 0) tZP1CorrHG[quad]  += (Float_t) (rawData.GetADCValue()-meanPed[pedindex]); 
234         else tZP1CorrLG[quad]  += (Float_t) (rawData.GetADCValue()-meanPed[pedindex+5]); 
235       }
236       else if(det == 3){ 
237         if(quad==1){     
238           pedindex = quad+20;
239           if(gain == 0) dZEMCorrHG += (Float_t) (rawData.GetADCValue()-meanPed[pedindex]); 
240           else dZEMCorrLG += (Float_t) (rawData.GetADCValue()-meanPed[pedindex+2]); 
241         }
242         else if(quad==2){ 
243           pedindex = rawData.GetSector(1)+21;
244           if(gain == 0) dZEMCorrHG += (Float_t) (rawData.GetADCValue()-meanPed[pedindex]); 
245           else dZEMCorrLG += (Float_t) (rawData.GetADCValue()-meanPed[pedindex+2]); 
246         }
247       }
248       else if(det == 4){       
249         pedindex = rawData.GetSector(1)+24;
250         if(gain == 0) tZN2CorrHG[quad]  += (Float_t) (rawData.GetADCValue()-meanPed[pedindex]); 
251         else tZN2CorrLG[quad]  += (Float_t) (rawData.GetADCValue()-meanPed[pedindex+2]); 
252       }
253       else if(det == 5){
254         pedindex = rawData.GetSector(1)+34;
255         if(gain == 0) tZP2CorrHG[quad]  += (Float_t) (rawData.GetADCValue()-meanPed[pedindex]); 
256         else tZP2CorrLG[quad]  += (Float_t) (rawData.GetADCValue()-meanPed[pedindex+5]); 
257       }
258     }
259   }
260     
261   // reconstruct the event
262     ReconstructEvent(clustersTree, tZN1CorrHG, tZP1CorrHG, tZN2CorrHG, 
263         tZP2CorrHG, tZN1CorrLG, tZP1CorrLG, tZN2CorrLG, 
264         tZP2CorrLG, dZEMCorrHG);
265
266 }
267
268 //_____________________________________________________________________________
269 void AliZDCReconstructor::ReconstructEvent(TTree *clustersTree, 
270                 Float_t* ZN1ADCCorrHG, Float_t* ZP1ADCCorrHG, 
271                 Float_t* ZN2ADCCorrHG, Float_t* ZP2ADCCorrHG, 
272                 Float_t* ZN1ADCCorrLG, Float_t* ZP1ADCCorrLG, 
273                 Float_t* ZN2ADCCorrLG, Float_t* ZP2ADCCorrLG, 
274                 Float_t corrADCZEMHG) const
275 {
276   // ***** Reconstruct one event
277   
278   // *** RECONSTRUCTION FROM SIMULATED DATA
279   // It passes trhough the no. of phe which is known from simulations
280   //  ---      ADCchannel -> photoelectrons
281   // NB-> PM gain = 10^(5), ADC resolution = 6.4*10^(-7)
282   // Move to V965 (E.S.,15/09/04) NB-> PM gain = 10^(5), ADC resolution = 8*10^(-7)
283   //Float_t zn1phe, zp1phe, zemphe, zn2phe, zp2phe, convFactor = 0.08;
284   //zn1phe  = ZN1Corr/convFactor;
285   //zp1phe  = ZP1Corr/convFactor;
286   //zemphe = ZEMCorr/convFactor;
287   //zn2phe  = ZN2Corr/convFactor;
288   //zp2phe  = ZP2Corr/convFactor;
289   ////if AliDebug(1,Form("\n    znphe = %f, zpphe = %f, zemphe = %f\n",znphe, zpphe, zemphe);
290   //
291   ////  ---      Energy calibration
292   //// Conversion factors for hadronic ZDCs goes from phe yield to TRUE 
293   //// incident energy (conversion from GeV to TeV is included); while for EM 
294   //// calos conversion is from light yield to detected energy calculated by
295   //// GEANT NB -> ZN and ZP conversion factors are constant since incident
296   //// spectators have all the same energy, ZEM energy is obtained through a
297   //// fit over the whole range of incident particle energies 
298   //// (obtained with full HIJING simulations) 
299   //Float_t zn1energy, zp1energy, zemenergy, zdc1energy, zn2energy, zp2energy, zdc2energy;
300   //Float_t zn1phexTeV=329., zp1phexTeV=369., zn2phexTeV=329., zp2phexTeV=369.;
301   //zn1energy  = zn1phe/zn1phexTeV;
302   //zp1energy  = zp1phe/zp1phexTeV;
303   //zdc1energy = zn1energy+zp1energy;
304   //zn2energy  = zn2phe/zn2phexTeV;
305   //zp2energy  = zp2phe/zp2phexTeV;
306   //zdc2energy = zn2energy+zp2energy;
307   //zemenergy = -4.81+0.3238*zemphe;
308   //if(zemenergy<0) zemenergy=0;
309   ////  if AliDebug(1,Form("    znenergy = %f TeV, zpenergy = %f TeV, zdcenergy = %f GeV, "
310   ////                     "\n          zemenergy = %f TeV\n", znenergy, zpenergy, 
311   ////                     zdcenergy, zemenergy);
312   ////  if(zdcenergy==0)
313   ////    if AliDebug(1,Form("\n\n      ###     ATTENZIONE!!! -> ev# %d: znenergy = %f TeV, zpenergy = %f TeV, zdcenergy = %f GeV, "
314   ////                       " zemenergy = %f TeV\n\n", fMerger->EvNum(), znenergy, zpenergy, zdcenergy, zemenergy); 
315   
316   //
317   // *** RECONSTRUCTION FROM "REAL" DATA
318   //
319   // Retrieving calibration data
320   Float_t equalCoeffZN1[5], equalCoeffZP1[5], equalCoeffZN2[5], equalCoeffZP2[5];
321   for(Int_t ji=0; ji<5; ji++){
322      equalCoeffZN1[ji] = fCalibData->GetZN1EqualCoeff(ji);
323      equalCoeffZP1[ji] = fCalibData->GetZP1EqualCoeff(ji); 
324      equalCoeffZN2[ji] = fCalibData->GetZN2EqualCoeff(ji); 
325      equalCoeffZP2[ji] = fCalibData->GetZP2EqualCoeff(ji); 
326   }
327   //
328   Float_t calibEne[4];
329   for(Int_t ij=0; ij<4; ij++) calibEne[ij] = fCalibData->GetEnCalib(ij);
330   //
331   Float_t endPointZEM = fCalibData->GetZEMEndValue();
332   Float_t cutFractionZEM = fCalibData->GetZEMCutFraction();
333   Float_t dZEMSup = fCalibData->GetDZEMSup();
334   Float_t dZEMInf = fCalibData->GetDZEMInf();
335   //
336   Float_t cutValueZEM = endPointZEM*cutFractionZEM;
337   Float_t supValueZEM = cutValueZEM+(endPointZEM*dZEMSup);
338   Float_t infValueZEM = cutValueZEM-(endPointZEM*dZEMInf);
339   //
340   Float_t maxValEZN1 = fCalibData->GetEZN1MaxValue();
341   Float_t maxValEZP1 = fCalibData->GetEZP1MaxValue();
342   Float_t maxValEZDC1 = fCalibData->GetEZDC1MaxValue();
343   Float_t maxValEZN2 = fCalibData->GetEZN2MaxValue();
344   Float_t maxValEZP2 = fCalibData->GetEZP2MaxValue();
345   Float_t maxValEZDC2 = fCalibData->GetEZDC2MaxValue();
346   //
347   //printf("\n\t AliZDCReconstructor -> ZEMEndPoint %1.0f, ZEMCutValue %1.0f,"
348   //   " ZEMSupValue %1.0f, ZEMInfValue %1.0f\n",endPointZEM,cutValueZEM,supValueZEM,infValueZEM);
349   
350   // Equalization of detector responses
351   Float_t equalTowZN1HG[5], equalTowZN2HG[5], equalTowZP1HG[5], equalTowZP2HG[5];
352   Float_t equalTowZN1LG[5], equalTowZN2LG[5], equalTowZP1LG[5], equalTowZP2LG[5];
353   for(Int_t gi=0; gi<5; gi++){
354      equalTowZN1HG[gi] = ZN1ADCCorrHG[gi]*equalCoeffZN1[gi];
355      equalTowZP1HG[gi] = ZP1ADCCorrHG[gi]*equalCoeffZP1[gi];
356      equalTowZN2HG[gi] = ZN2ADCCorrHG[gi]*equalCoeffZN2[gi];
357      equalTowZP2HG[gi] = ZP2ADCCorrHG[gi]*equalCoeffZP2[gi];
358      //
359      equalTowZN1LG[gi] = ZN1ADCCorrLG[gi]*equalCoeffZN1[gi];
360      equalTowZP1LG[gi] = ZP1ADCCorrLG[gi]*equalCoeffZP1[gi];
361      equalTowZN2LG[gi] = ZN2ADCCorrLG[gi]*equalCoeffZN2[gi];
362      equalTowZP2LG[gi] = ZP2ADCCorrLG[gi]*equalCoeffZP2[gi];
363   }
364   
365   // Energy calibration of detector responses
366   Float_t calibTowZN1HG[5], calibTowZN2HG[5], calibTowZP1HG[5], calibTowZP2HG[5];
367   Float_t calibSumZN1HG=0., calibSumZN2HG=0., calibSumZP1HG=0., calibSumZP2HG=0.;
368   Float_t calibTowZN1LG[5], calibTowZN2LG[5], calibTowZP1LG[5], calibTowZP2LG[5];
369   Float_t calibSumZN1LG=0., calibSumZN2LG=0., calibSumZ12LG=0., calibSumZP2LG=0.;
370   for(Int_t gi=0; gi<5; gi++){
371      calibTowZN1HG[gi] = equalTowZN1HG[gi]*calibEne[0];
372      calibTowZP1HG[gi] = equalTowZP1HG[gi]*calibEne[1];
373      calibTowZN2HG[gi] = equalTowZN2HG[gi]*calibEne[2];
374      calibTowZP2HG[gi] = equalTowZP2HG[gi]*calibEne[3];
375      calibSumZN1HG += calibTowZN1HG[gi];
376      calibSumZP1HG += calibTowZP1HG[gi];
377      calibSumZN2HG += calibTowZN2HG[gi];
378      calibSumZP2HG += calibTowZP2HG[gi];
379      //
380      calibTowZN1LG[gi] = equalTowZN1LG[gi]*calibEne[0];
381      calibTowZP1LG[gi] = equalTowZP1LG[gi]*calibEne[1];
382      calibTowZN2LG[gi] = equalTowZN2LG[gi]*calibEne[2];
383      calibTowZP2LG[gi] = equalTowZP2LG[gi]*calibEne[3];
384      calibSumZN1LG += calibTowZN1LG[gi];
385      calibSumZ12LG += calibTowZP1LG[gi];
386      calibSumZN2LG += calibTowZN2LG[gi];
387      calibSumZP2LG += calibTowZP2LG[gi];
388   }
389   
390   //  ---      Number of detected spectator nucleons
391   //  *** N.B. -> It works only in Pb-Pb
392   Int_t nDetSpecNLeft, nDetSpecPLeft, nDetSpecNRight, nDetSpecPRight;
393   nDetSpecNLeft = (Int_t) (calibSumZN1HG/2.760);
394   nDetSpecPLeft = (Int_t) (calibSumZP1HG/2.760);
395   nDetSpecNRight = (Int_t) (calibSumZN2HG/2.760);
396   nDetSpecPRight = (Int_t) (calibSumZP2HG/2.760);
397   /*printf("\n\t AliZDCReconstructor -> nDetSpecNLeft %d, nDetSpecPLeft %d,"
398     " nDetSpecNRight %d, nDetSpecPRight %d\n",nDetSpecNLeft, nDetSpecPLeft, 
399     nDetSpecNRight, nDetSpecPRight);*/
400
401   //  ---      Number of generated spectator nucleons (from HIJING parameterization)
402   Int_t nGenSpecNLeft=0, nGenSpecPLeft=0, nGenSpecLeft=0;
403   Int_t nGenSpecNRight=0, nGenSpecPRight=0, nGenSpecRight=0;
404   Double_t impPar=0.;
405   //
406   // *** RECONSTRUCTION FROM SIMULATED DATA
407   // Cut value for Ezem (GeV)
408   // ### Results from production  -> 0<b<18 fm (Apr 2002)
409   /*Float_t eZEMCut = 420.;
410   Float_t deltaEZEMSup = 690.; 
411   Float_t deltaEZEMInf = 270.; 
412   if(zemenergy > (eZEMCut+deltaEZEMSup)){
413     nGenSpecNLeft  = (Int_t) (fZNCen->Eval(ZN1CalibSum));
414     nGenSpecPLeft  = (Int_t) (fZPCen->Eval(ZP1CalibSum));
415     nGenSpecLeft   = (Int_t) (fZDCCen->Eval(ZN1CalibSum+ZP1CalibSum));
416     nGenSpecNRight = (Int_t) (fZNCen->Eval(ZN2CalibSum));
417     nGenSpecPRight = (Int_t) (fZNCen->Eval(ZP2CalibSum));
418     nGenSpecRight  = (Int_t) (fZNCen->Eval(ZN2CalibSum+ZP2CalibSum));
419     impPar  = fbCen->Eval(ZN1CalibSum+ZP1CalibSum);
420   }
421   else if(zemenergy < (eZEMCut-deltaEZEMInf)){
422     nGenSpecNLeft = (Int_t) (fZNPer->Eval(ZN1CalibSum)); 
423     nGenSpecPLeft = (Int_t) (fZPPer->Eval(ZP1CalibSum));
424     nGenSpecLeft  = (Int_t) (fZDCPer->Eval(ZN1CalibSum+ZP1CalibSum));
425     impPar   = fbPer->Eval(ZN1CalibSum+ZP1CalibSum);
426   }
427   else if(zemenergy >= (eZEMCut-deltaEZEMInf) && zemenergy <= (eZEMCut+deltaEZEMSup)){
428     nGenSpecNLeft = (Int_t) (fZEMn->Eval(zemenergy));
429     nGenSpecPLeft = (Int_t) (fZEMp->Eval(zemenergy));
430     nGenSpecLeft  = (Int_t)(fZEMsp->Eval(zemenergy));
431     impPar   =  fZEMb->Eval(zemenergy);
432   }
433   // ### Results from production  -> 0<b<18 fm (Apr 2002)
434   if(ZN1CalibSum>162.)  nGenSpecNLeft = (Int_t) (fZEMn->Eval(zemenergy));
435   if(ZP1CalibSum>59.75)  nGenSpecPLeft = (Int_t) (fZEMp->Eval(zemenergy));
436   if(ZN1CalibSum+ZP1CalibSum>221.5) nGenSpecLeft  = (Int_t)(fZEMsp->Eval(zemenergy));
437   if(ZN1CalibSum+ZP1CalibSum>220.)  impPar    =  fZEMb->Eval(zemenergy);
438   */
439   //
440   //
441   // *** RECONSTRUCTION FROM REAL DATA
442   //
443   if(corrADCZEMHG > supValueZEM){
444     nGenSpecNLeft  = (Int_t) (fZNCen->Eval(calibSumZN1HG));
445     nGenSpecPLeft  = (Int_t) (fZPCen->Eval(calibSumZP1HG));
446     nGenSpecLeft   = (Int_t) (fZDCCen->Eval(calibSumZN1HG+calibSumZP1HG));
447     nGenSpecNRight = (Int_t) (fZNCen->Eval(calibSumZN2HG));
448     nGenSpecPRight = (Int_t) (fZNCen->Eval(calibSumZP2HG));
449     nGenSpecRight  = (Int_t) (fZNCen->Eval(calibSumZN2HG+calibSumZP2HG));
450     impPar  = fbCen->Eval(calibSumZN1HG+calibSumZP1HG);
451   }
452   else if(corrADCZEMHG < infValueZEM){
453     nGenSpecNLeft = (Int_t) (fZNPer->Eval(calibSumZN1HG)); 
454     nGenSpecPLeft = (Int_t) (fZPPer->Eval(calibSumZP1HG));
455     nGenSpecLeft  = (Int_t) (fZDCPer->Eval(calibSumZN1HG+calibSumZP1HG));
456     impPar   = fbPer->Eval(calibSumZN1HG+calibSumZP1HG);
457   }
458   else if(corrADCZEMHG >= infValueZEM && corrADCZEMHG <= supValueZEM){
459     nGenSpecNLeft = (Int_t) (fZEMn->Eval(corrADCZEMHG));
460     nGenSpecPLeft = (Int_t) (fZEMp->Eval(corrADCZEMHG));
461     nGenSpecLeft  = (Int_t)(fZEMsp->Eval(corrADCZEMHG));
462     impPar   =  fZEMb->Eval(corrADCZEMHG);
463   }
464   // 
465   if(calibSumZN1HG/maxValEZN1>1.)  nGenSpecNLeft = (Int_t) (fZEMn->Eval(corrADCZEMHG));
466   if(calibSumZP1HG/maxValEZP1>1.)  nGenSpecPLeft = (Int_t) (fZEMp->Eval(corrADCZEMHG));
467   if((calibSumZN1HG+calibSumZP1HG/maxValEZDC1)>1.){
468      nGenSpecLeft = (Int_t)(fZEMsp->Eval(corrADCZEMHG));
469      impPar = fZEMb->Eval(corrADCZEMHG);
470   }
471   if(calibSumZN2HG/maxValEZN2>1.)  nGenSpecNRight = (Int_t) (fZEMn->Eval(corrADCZEMHG));
472   if(calibSumZP2HG/maxValEZP2>1.)  nGenSpecPRight = (Int_t) (fZEMp->Eval(corrADCZEMHG));
473   if((calibSumZN2HG+calibSumZP2HG/maxValEZDC2)>1.) nGenSpecRight = (Int_t)(fZEMsp->Eval(corrADCZEMHG));
474   //
475   if(nGenSpecNLeft>125)    nGenSpecNLeft=125;
476   else if(nGenSpecNLeft<0) nGenSpecNLeft=0;
477   if(nGenSpecPLeft>82)     nGenSpecPLeft=82;
478   else if(nGenSpecPLeft<0) nGenSpecPLeft=0;
479   if(nGenSpecLeft>207)     nGenSpecLeft=207;
480   else if(nGenSpecLeft<0)  nGenSpecLeft=0;
481   
482   //  ---      Number of generated participants (from HIJING parameterization)
483   Int_t nPart, nPartTotLeft, nPartTotRight;
484   nPart = 207-nGenSpecNLeft-nGenSpecPLeft;
485   nPartTotLeft = 207-nGenSpecLeft;
486   nPartTotRight = 207-nGenSpecRight;
487   if(nPart<0) nPart=0;
488   if(nPartTotLeft<0) nPartTotLeft=0;
489   if(nPartTotRight<0) nPartTotRight=0;
490   //
491   // *** DEBUG ***
492   printf("\n\t AliZDCReconstructor -> calibSumZN1HG %1.0f, calibSumZP1HG %1.0f,"
493       "  calibSumZN2HG %1.0f, calibSumZP2HG %1.0f, corrADCZEMHG %1.0f\n", 
494       calibSumZN1HG,calibSumZP1HG,calibSumZN2HG,calibSumZP2HG,corrADCZEMHG);
495   printf("\t AliZDCReconstructor -> nGenSpecNLeft %d, nGenSpecPLeft %d, nGenSpecLeft %d\n"
496       "\t\t nGenSpecNRight %d, nGenSpecPRight %d, nGenSpecRight %d\n", 
497       nGenSpecNLeft, nGenSpecPLeft, nGenSpecLeft, 
498       nGenSpecNRight, nGenSpecPRight, nGenSpecRight);
499   printf("\t AliZDCReconstructor ->  NpartL %d,  NpartR %d,  b %1.2f fm\n\n",nPartTotLeft, nPartTotRight, impPar);
500
501   // create the output tree
502   AliZDCReco reco(calibSumZN1HG, calibSumZP1HG, calibSumZN2HG, calibSumZP2HG, 
503                   calibTowZN1LG, calibTowZN2LG, calibTowZP1LG, calibTowZP2LG, 
504                   corrADCZEMHG, 
505                   nDetSpecNLeft, nDetSpecPLeft, nDetSpecNRight, nDetSpecPRight, 
506                   nGenSpecNLeft, nGenSpecPLeft, nGenSpecLeft, nGenSpecNRight, 
507                   nGenSpecPRight, nGenSpecRight,
508                   nPartTotLeft, nPartTotRight, impPar);
509                   
510   AliZDCReco* preco = &reco;
511   const Int_t kBufferSize = 4000;
512   clustersTree->Branch("ZDC", "AliZDCReco", &preco, kBufferSize);
513
514   // write the output tree
515   clustersTree->Fill();
516 }
517
518 //_____________________________________________________________________________
519 void AliZDCReconstructor::FillZDCintoESD(TTree *clustersTree, AliESDEvent* esd) const
520 {
521   // fill energies and number of participants to the ESD
522
523   AliZDCReco reco;
524   AliZDCReco* preco = &reco;
525   clustersTree->SetBranchAddress("ZDC", &preco);
526
527   clustersTree->GetEntry(0);
528   /*Double_t tZN1Ene[4], tZN2Ene[4];
529   for(Int_t i=0; i<4; i++){
530      tZN1Ene[i] = reco.GetZN1EnTow(i);
531      tZN2Ene[i] = reco.GetZN2EnTow(i);
532   }
533   esd->SetZDC(tZN1Ene, tZN2Ene, reco.GetZN1Energy(), reco.GetZP1Energy(), reco.GetZEMsignal(),
534               reco.GetZN2Energy(), reco.GetZP2Energy(), 
535               reco.GetNPartLeft());
536   */
537   esd->SetZDC(reco.GetZN1Energy(), reco.GetZP1Energy(), reco.GetZEMsignal(),
538               reco.GetZN2Energy(), reco.GetZP2Energy(), 
539               reco.GetNPartLeft());
540   
541   /*Double_t tZN1Ene[4], tZN2Ene[4];
542   for(Int_t i=0; i<4; i++){
543      tZN1Ene[i] = reco.GetZN1EnTow(i);
544      tZN2Ene[i] = reco.GetZN2EnTow(i);
545   }*/
546   //esd->SetZN1TowerEnergy(tZN1Ene);
547   //esd->SetZN2TowerEnergy(tZN2Ene);
548   
549 }
550
551 //_____________________________________________________________________________
552 AliCDBStorage* AliZDCReconstructor::SetStorage(const char *uri) 
553 {
554   // Setting the storage
555
556   Bool_t deleteManager = kFALSE;
557   
558   AliCDBManager *manager = AliCDBManager::Instance();
559   AliCDBStorage *defstorage = manager->GetDefaultStorage();
560   
561   if(!defstorage || !(defstorage->Contains("ZDC"))){ 
562      AliWarning("No default storage set or default storage doesn't contain ZDC!");
563      manager->SetDefaultStorage(uri);
564      deleteManager = kTRUE;
565   }
566  
567   AliCDBStorage *storage = manager->GetDefaultStorage();
568
569   if(deleteManager){
570     AliCDBManager::Instance()->UnsetDefaultStorage();
571     defstorage = 0;   // the storage is killed by AliCDBManager::Instance()->Destroy()
572   }
573
574   return storage; 
575 }
576
577 //_____________________________________________________________________________
578 AliZDCCalibData* AliZDCReconstructor::GetCalibData() const
579 {
580
581   // Getting calibration object for ZDC set
582
583   AliCDBEntry  *entry = AliCDBManager::Instance()->Get("ZDC/Calib/Data");
584   if(!entry) AliFatal("No calibration data loaded!");  
585
586   AliZDCCalibData *calibdata = dynamic_cast<AliZDCCalibData*>  (entry->GetObject());
587   if(!calibdata)  AliFatal("Wrong calibration object in calibration  file!");
588
589   return calibdata;
590 }
591 /**************************************************************************
592  * Copyright(c) 1998-1999, ALICE Experiment at CERN, All rights reserved. *
593  *                                                                        *
594  * Author: The ALICE Off-line Project.                                    *
595  * Contributors are mentioned in the code where appropriate.              *
596  *                                                                        *
597  * Permission to use, copy, modify and distribute this software and its   *
598  * documentation strictly for non-commercial purposes is hereby granted   *
599  * without fee, provided that the above copyright notice appears in all   *
600  * copies and that both the copyright notice and this permission notice   *
601  * appear in the supporting documentation. The authors make no claims     *
602  * about the suitability of this software for any purpose. It is          *
603  * provided "as is" without express or implied warranty.                  *
604  **************************************************************************/
605
606 /* $Id$ */
607
608 ///////////////////////////////////////////////////////////////////////////////
609 //                                                                           //
610 // class for ZDC reconstruction                                              //
611 //                                                                           //
612 ///////////////////////////////////////////////////////////////////////////////
613
614
615 #include <TF1.h>
616
617 #include "AliRunLoader.h"
618 #include "AliRawReader.h"
619 #include "AliESDEvent.h"
620 #include "AliZDCDigit.h"
621 #include "AliZDCRawStream.h"
622 #include "AliZDCReco.h"
623 #include "AliZDCReconstructor.h"
624 #include "AliZDCCalibData.h"
625
626
627 ClassImp(AliZDCReconstructor)
628
629
630 //_____________________________________________________________________________
631 AliZDCReconstructor:: AliZDCReconstructor() :
632
633   fZNCen(new TF1("fZNCen", 
634         "(-2.287920+sqrt(2.287920*2.287920-4*(-0.007629)*(11.921710-x)))/(2*(-0.007629))",0.,164.)),
635   fZNPer(new TF1("fZNPer",
636       "(-37.812280-sqrt(37.812280*37.812280-4*(-0.190932)*(-1709.249672-x)))/(2*(-0.190932))",0.,164.)),
637   fZPCen(new TF1("fZPCen",
638        "(-1.321353+sqrt(1.321353*1.321353-4*(-0.007283)*(3.550697-x)))/(2*(-0.007283))",0.,60.)),
639   fZPPer(new TF1("fZPPer",
640       "(-42.643308-sqrt(42.643308*42.643308-4*(-0.310786)*(-1402.945615-x)))/(2*(-0.310786))",0.,60.)),
641   fZDCCen(new TF1("fZDCCen",
642       "(-1.934991+sqrt(1.934991*1.934991-4*(-0.004080)*(15.111124-x)))/(2*(-0.004080))",0.,225.)),
643   fZDCPer(new TF1("fZDCPer",
644       "(-34.380639-sqrt(34.380639*34.380639-4*(-0.104251)*(-2612.189017-x)))/(2*(-0.104251))",0.,225.)),
645   fbCen(new TF1("fbCen","-0.056923+0.079703*x-0.0004301*x*x+0.000001366*x*x*x",0.,220.)),
646   fbPer(new TF1("fbPer","17.943998-0.046846*x+0.000074*x*x",0.,220.)),
647   //
648   fZEMn(new TF1("fZEMn","121.7-0.1934*x+0.00007565*x*x",0.,1200.)),
649   fZEMp(new TF1("fZEMp","80.05-0.1315*x+0.00005327*x*x",0.,1200.)),
650   fZEMsp(new TF1("fZEMsp","201.7-0.325*x+0.0001292*x*x",0.,1200.)),
651   fZEMb(new TF1("fZEMb",
652         "13.83-0.02851*x+5.101e-5*x*x-7.305e-8*x*x*x+5.101e-11*x*x*x*x-1.25e-14*x*x*x*x*x",0.,1200.)),
653   //
654   fCalibData(GetCalibData())
655
656 {
657   // **** Default constructor
658
659 }
660
661
662 //_____________________________________________________________________________
663 AliZDCReconstructor::~AliZDCReconstructor()
664 {
665 // destructor
666
667   delete fZNCen;
668   delete fZNPer;
669   delete fZPCen;
670   delete fZPPer;
671   delete fZDCCen;
672   delete fZDCPer;
673   delete fbCen;
674   delete fbPer;
675   delete fZEMn;
676   delete fZEMp;
677   delete fZEMsp;
678   delete fZEMb;
679
680 }
681
682
683 //_____________________________________________________________________________
684 void AliZDCReconstructor::Reconstruct(TTree* digitsTree, TTree* clustersTree) const
685 {
686   // *** Local ZDC reconstruction for digits
687   // Works on the current event
688     
689   // Retrieving calibration data  
690   Float_t meanPed[47];
691   for(Int_t jj=0; jj<47; jj++) meanPed[jj] = fCalibData->GetMeanPed(jj);
692
693   // get digits
694   AliZDCDigit digit;
695   AliZDCDigit* pdigit = &digit;
696   digitsTree->SetBranchAddress("ZDC", &pdigit);
697
698   // loop over digits
699   Float_t tZN1CorrHG[]={0.,0.,0.,0.,0.}, tZP1CorrHG[]={0.,0.,0.,0.,0.}; 
700   Float_t dZEMCorrHG=0.; 
701   Float_t tZN2CorrHG[]={0.,0.,0.,0.,0.}, tZP2CorrHG[]={0.,0.,0.,0.,0.};
702   Float_t tZN1CorrLG[]={0.,0.,0.,0.,0.}, tZP1CorrLG[]={0.,0.,0.,0.,0.};
703   Float_t dZEMCorrLG=0.; 
704   Float_t tZN2CorrLG[]={0.,0.,0.,0.,0.}, tZP2CorrLG[]={0.,0.,0.,0.,0.};
705   
706   //printf("\n\t # of digits in tree: %d\n",(Int_t) digitsTree->GetEntries());
707   for (Int_t iDigit = 0; iDigit < (digitsTree->GetEntries()/2); iDigit++) {
708     digitsTree->GetEntry(iDigit);
709     if (!pdigit) continue;
710     //pdigit->Print("");
711     //  
712     Int_t det = digit.GetSector(0);
713     Int_t quad = digit.GetSector(1);
714     Int_t pedindex = -1;
715     //printf("\n\t #%d det %d quad %d", iDigit, det, quad);
716     //
717     if(det == 1){ // *** ZN1
718        pedindex = quad;
719        tZN1CorrHG[quad] = (Float_t) (digit.GetADCValue(0)-meanPed[pedindex]);
720        if(tZN1CorrHG[quad]<0.) tZN1CorrHG[quad] = 0.;
721        tZN1CorrLG[quad] = (Float_t) (digit.GetADCValue(1)-meanPed[pedindex+5]);
722        if(tZN1CorrLG[quad]<0.) tZN1CorrLG[quad] = 0.;
723        //printf("\t pedindex %d tZN1CorrHG[%d] = %1.0f tZN1CorrLG[%d] = %1.0f", 
724        //       pedindex, quad, tZN1CorrHG[quad], quad, tZN1CorrLG[quad]);
725     }
726     else if(det == 2){ // *** ZP1
727        pedindex = quad+10;
728        tZP1CorrHG[quad] = (Float_t) (digit.GetADCValue(0)-meanPed[pedindex]);
729        if(tZP1CorrLG[quad]<0.) tZP1CorrLG[quad] = 0.;
730        tZP1CorrLG[quad] = (Float_t) (digit.GetADCValue(1)-meanPed[pedindex+5]);
731        if(tZP1CorrHG[quad]<0.) tZP1CorrHG[quad] = 0.;
732        //printf("\t pedindex %d tZP1CorrHG[%d] = %1.0f tZP1CorrLG[%d] = %1.0f", 
733        //       pedindex, quad, tZP1CorrHG[quad], quad, tZP1CorrLG[quad]);
734     }
735     else if(det == 3){
736        if(quad == 1){       // *** ZEM1  
737          pedindex = quad+19;
738          dZEMCorrHG += (Float_t) (digit.GetADCValue(0)-meanPed[pedindex]); 
739          if(dZEMCorrHG<0.) dZEMCorrHG = 0.;
740          dZEMCorrLG += (Float_t) (digit.GetADCValue(1)-meanPed[pedindex+2]); 
741          if(dZEMCorrLG<0.) dZEMCorrLG = 0.;
742          //printf("\t pedindex %d ADC(0) = %d ped = %1.0f ADCCorr = %1.0f\n", 
743          //     pedindex, digit.GetADCValue(0), meanPed[pedindex], dZEMCorrHG);
744          //printf("\t pedindex %d ADC(1) = %d ped = %1.0f ADCCorr = %1.0f\n", 
745          //     pedindex+2, digit.GetADCValue(1), meanPed[pedindex+2], dZEMCorrLG);
746          ////printf("\t pedindex %d dZEMCorrHG = %1.0f dZEMCorrLG = %1.0f\n", pedindex, dZEMCorrHG, dZEMCorrLG);
747        }
748        else if(quad == 2){  // *** ZEM1
749          pedindex = quad+19;
750          dZEMCorrHG += (Float_t) (digit.GetADCValue(0)-meanPed[pedindex]); 
751          if(dZEMCorrHG<0.) dZEMCorrHG = 0.;
752          dZEMCorrLG += (Float_t) (digit.GetADCValue(1)-meanPed[pedindex+2]); 
753          if(dZEMCorrLG<0.) dZEMCorrLG = 0.;
754          //printf("\t pedindex %d ADC(0) = %d ped = %1.0f ADCCorr = %1.0f\n", 
755          //     pedindex, digit.GetADCValue(0), meanPed[pedindex], dZEMCorrHG);
756          //printf("\t pedindex %d ADC(1) = %d ped = %1.0f ADCCorr = %1.0f\n", 
757          //     pedindex+2, digit.GetADCValue(1),meanPed[pedindex+2], dZEMCorrLG);
758          ////printf("\t pedindex %d dZEMCorrHG = %1.0f dZEMCorrLG = %1.0f\n", pedindex, dZEMCorrHG, dZEMCorrLG);
759        }
760     }
761     else if(det == 4){  // *** ZN2
762        pedindex = quad+24;
763        tZN2CorrHG[quad] = (Float_t) (digit.GetADCValue(0)-meanPed[pedindex]);
764        if(tZN2CorrHG[quad]<0.) tZN2CorrHG[quad] = 0.;
765        tZN2CorrLG[quad] = (Float_t) (digit.GetADCValue(1)-meanPed[pedindex+5]);
766        if(tZN2CorrLG[quad]<0.) tZN2CorrLG[quad] = 0.;
767        //printf("\t pedindex %d tZN2CorrHG[%d] = %1.0f tZN2CorrLG[%d] = %1.0f\n", 
768        //       pedindex, quad, tZN2CorrHG[quad], quad, tZN2CorrLG[quad]);
769     }
770     else if(det == 5){  // *** ZP2 
771        pedindex = quad+34;
772        tZP2CorrHG[quad] = (Float_t) (digit.GetADCValue(0)-meanPed[pedindex]);
773        if(tZP2CorrHG[quad]<0.) tZP2CorrHG[quad] = 0.;
774        tZP2CorrLG[quad] = (Float_t) (digit.GetADCValue(1)-meanPed[pedindex+5]);
775        if(tZP2CorrLG[quad]<0.) tZP2CorrLG[quad] = 0.;
776        //printf("\t pedindex %d tZP2CorrHG[%d] = %1.0f tZP2CorrLG[%d] = %1.0f\n", 
777        //       pedindex, quad, tZP2CorrHG[quad], quad, tZP2CorrLG[quad]);
778     }
779   }
780
781   // reconstruct the event
782     ReconstructEvent(clustersTree, tZN1CorrHG, tZP1CorrHG, tZN2CorrHG, 
783         tZP2CorrHG, tZN1CorrLG, tZP1CorrLG, tZN2CorrLG, 
784         tZP2CorrLG, dZEMCorrHG);
785
786 }
787
788 //_____________________________________________________________________________
789 void AliZDCReconstructor::Reconstruct(AliRawReader* rawReader, TTree* clustersTree) const
790 {
791   // *** ZDC raw data reconstruction
792   // Works on the current event
793   
794   // Retrieving calibration data  
795   Float_t meanPed[47];
796   for(Int_t jj=0; jj<47; jj++) meanPed[jj] = fCalibData->GetMeanPed(jj);
797
798   rawReader->Reset();
799
800   // loop over raw data rawDatas
801   Float_t tZN1CorrHG[]={0.,0.,0.,0.,0.}, tZP1CorrHG[]={0.,0.,0.,0.,0.};
802   Float_t dZEMCorrHG=0.;
803   Float_t tZN2CorrHG[]={0.,0.,0.,0.,0.}, tZP2CorrHG[]={0.,0.,0.,0.,0.};
804   Float_t tZN1CorrLG[]={0.,0.,0.,0.,0.}, tZP1CorrLG[]={0.,0.,0.,0.,0.};
805   Float_t dZEMCorrLG=0.; 
806   Float_t tZN2CorrLG[]={0.,0.,0.,0.,0.}, tZP2CorrLG[]={0.,0.,0.,0.,0.};
807   //
808   AliZDCRawStream rawData(rawReader);
809   while (rawData.Next()) {
810     if(rawData.IsADCDataWord()){
811       Int_t det = rawData.GetSector(0);
812       Int_t quad = rawData.GetSector(1);
813       Int_t gain = rawData.GetADCGain();
814       Int_t pedindex;
815       //
816       if(det == 1){    
817         pedindex = quad;
818         if(gain == 0) tZN1CorrHG[quad]  += (Float_t) (rawData.GetADCValue()-meanPed[pedindex]); 
819         else tZN1CorrLG[quad]  += (Float_t) (rawData.GetADCValue()-meanPed[pedindex+5]); 
820       }
821       else if(det == 2){ 
822         pedindex = quad+10;
823         if(gain == 0) tZP1CorrHG[quad]  += (Float_t) (rawData.GetADCValue()-meanPed[pedindex]); 
824         else tZP1CorrLG[quad]  += (Float_t) (rawData.GetADCValue()-meanPed[pedindex+5]); 
825       }
826       else if(det == 3){ 
827         if(quad==1){     
828           pedindex = quad+20;
829           if(gain == 0) dZEMCorrHG += (Float_t) (rawData.GetADCValue()-meanPed[pedindex]); 
830           else dZEMCorrLG += (Float_t) (rawData.GetADCValue()-meanPed[pedindex+2]); 
831         }
832         else if(quad==2){ 
833           pedindex = rawData.GetSector(1)+21;
834           if(gain == 0) dZEMCorrHG += (Float_t) (rawData.GetADCValue()-meanPed[pedindex]); 
835           else dZEMCorrLG += (Float_t) (rawData.GetADCValue()-meanPed[pedindex+2]); 
836         }
837       }
838       else if(det == 4){       
839         pedindex = rawData.GetSector(1)+24;
840         if(gain == 0) tZN2CorrHG[quad]  += (Float_t) (rawData.GetADCValue()-meanPed[pedindex]); 
841         else tZN2CorrLG[quad]  += (Float_t) (rawData.GetADCValue()-meanPed[pedindex+2]); 
842       }
843       else if(det == 5){
844         pedindex = rawData.GetSector(1)+34;
845         if(gain == 0) tZP2CorrHG[quad]  += (Float_t) (rawData.GetADCValue()-meanPed[pedindex]); 
846         else tZP2CorrLG[quad]  += (Float_t) (rawData.GetADCValue()-meanPed[pedindex+5]); 
847       }
848     }
849   }
850     
851   // reconstruct the event
852     ReconstructEvent(clustersTree, tZN1CorrHG, tZP1CorrHG, tZN2CorrHG, 
853         tZP2CorrHG, tZN1CorrLG, tZP1CorrLG, tZN2CorrLG, 
854         tZP2CorrLG, dZEMCorrHG);
855
856 }
857
858 //_____________________________________________________________________________
859 void AliZDCReconstructor::ReconstructEvent(TTree *clustersTree, 
860                 Float_t* ZN1ADCCorrHG, Float_t* ZP1ADCCorrHG, 
861                 Float_t* ZN2ADCCorrHG, Float_t* ZP2ADCCorrHG, 
862                 Float_t* ZN1ADCCorrLG, Float_t* ZP1ADCCorrLG, 
863                 Float_t* ZN2ADCCorrLG, Float_t* ZP2ADCCorrLG, 
864                 Float_t corrADCZEMHG) const
865 {
866   // ***** Reconstruct one event
867   
868   // *** RECONSTRUCTION FROM SIMULATED DATA
869   // It passes trhough the no. of phe which is known from simulations
870   //  ---      ADCchannel -> photoelectrons
871   // NB-> PM gain = 10^(5), ADC resolution = 6.4*10^(-7)
872   // Move to V965 (E.S.,15/09/04) NB-> PM gain = 10^(5), ADC resolution = 8*10^(-7)
873   //Float_t zn1phe, zp1phe, zemphe, zn2phe, zp2phe, convFactor = 0.08;
874   //zn1phe  = ZN1Corr/convFactor;
875   //zp1phe  = ZP1Corr/convFactor;
876   //zemphe = ZEMCorr/convFactor;
877   //zn2phe  = ZN2Corr/convFactor;
878   //zp2phe  = ZP2Corr/convFactor;
879   ////if AliDebug(1,Form("\n    znphe = %f, zpphe = %f, zemphe = %f\n",znphe, zpphe, zemphe);
880   //
881   ////  ---      Energy calibration
882   //// Conversion factors for hadronic ZDCs goes from phe yield to TRUE 
883   //// incident energy (conversion from GeV to TeV is included); while for EM 
884   //// calos conversion is from light yield to detected energy calculated by
885   //// GEANT NB -> ZN and ZP conversion factors are constant since incident
886   //// spectators have all the same energy, ZEM energy is obtained through a
887   //// fit over the whole range of incident particle energies 
888   //// (obtained with full HIJING simulations) 
889   //Float_t zn1energy, zp1energy, zemenergy, zdc1energy, zn2energy, zp2energy, zdc2energy;
890   //Float_t zn1phexTeV=329., zp1phexTeV=369., zn2phexTeV=329., zp2phexTeV=369.;
891   //zn1energy  = zn1phe/zn1phexTeV;
892   //zp1energy  = zp1phe/zp1phexTeV;
893   //zdc1energy = zn1energy+zp1energy;
894   //zn2energy  = zn2phe/zn2phexTeV;
895   //zp2energy  = zp2phe/zp2phexTeV;
896   //zdc2energy = zn2energy+zp2energy;
897   //zemenergy = -4.81+0.3238*zemphe;
898   //if(zemenergy<0) zemenergy=0;
899   ////  if AliDebug(1,Form("    znenergy = %f TeV, zpenergy = %f TeV, zdcenergy = %f GeV, "
900   ////                     "\n          zemenergy = %f TeV\n", znenergy, zpenergy, 
901   ////                     zdcenergy, zemenergy);
902   ////  if(zdcenergy==0)
903   ////    if AliDebug(1,Form("\n\n      ###     ATTENZIONE!!! -> ev# %d: znenergy = %f TeV, zpenergy = %f TeV, zdcenergy = %f GeV, "
904   ////                       " zemenergy = %f TeV\n\n", fMerger->EvNum(), znenergy, zpenergy, zdcenergy, zemenergy); 
905   
906   //
907   // *** RECONSTRUCTION FROM "REAL" DATA
908   //
909   // Retrieving calibration data
910   Float_t equalCoeffZN1[5], equalCoeffZP1[5], equalCoeffZN2[5], equalCoeffZP2[5];
911   for(Int_t ji=0; ji<5; ji++){
912      equalCoeffZN1[ji] = fCalibData->GetZN1EqualCoeff(ji);
913      equalCoeffZP1[ji] = fCalibData->GetZP1EqualCoeff(ji); 
914      equalCoeffZN2[ji] = fCalibData->GetZN2EqualCoeff(ji); 
915      equalCoeffZP2[ji] = fCalibData->GetZP2EqualCoeff(ji); 
916   }
917   //
918   Float_t calibEne[4];
919   for(Int_t ij=0; ij<4; ij++) calibEne[ij] = fCalibData->GetEnCalib(ij);
920   //
921   Float_t endPointZEM = fCalibData->GetZEMEndValue();
922   Float_t cutFractionZEM = fCalibData->GetZEMCutFraction();
923   Float_t dZEMSup = fCalibData->GetDZEMSup();
924   Float_t dZEMInf = fCalibData->GetDZEMInf();
925   //
926   Float_t cutValueZEM = endPointZEM*cutFractionZEM;
927   Float_t supValueZEM = cutValueZEM+(endPointZEM*dZEMSup);
928   Float_t infValueZEM = cutValueZEM-(endPointZEM*dZEMInf);
929   //
930   Float_t maxValEZN1 = fCalibData->GetEZN1MaxValue();
931   Float_t maxValEZP1 = fCalibData->GetEZP1MaxValue();
932   Float_t maxValEZDC1 = fCalibData->GetEZDC1MaxValue();
933   Float_t maxValEZN2 = fCalibData->GetEZN2MaxValue();
934   Float_t maxValEZP2 = fCalibData->GetEZP2MaxValue();
935   Float_t maxValEZDC2 = fCalibData->GetEZDC2MaxValue();
936   //
937   //printf("\n\t AliZDCReconstructor -> ZEMEndPoint %1.0f, ZEMCutValue %1.0f,"
938   //   " ZEMSupValue %1.0f, ZEMInfValue %1.0f\n",endPointZEM,cutValueZEM,supValueZEM,infValueZEM);
939   
940   // Equalization of detector responses
941   Float_t equalTowZN1HG[5], equalTowZN2HG[5], equalTowZP1HG[5], equalTowZP2HG[5];
942   Float_t equalTowZN1LG[5], equalTowZN2LG[5], equalTowZP1LG[5], equalTowZP2LG[5];
943   for(Int_t gi=0; gi<5; gi++){
944      equalTowZN1HG[gi] = ZN1ADCCorrHG[gi]*equalCoeffZN1[gi];
945      equalTowZP1HG[gi] = ZP1ADCCorrHG[gi]*equalCoeffZP1[gi];
946      equalTowZN2HG[gi] = ZN2ADCCorrHG[gi]*equalCoeffZN2[gi];
947      equalTowZP2HG[gi] = ZP2ADCCorrHG[gi]*equalCoeffZP2[gi];
948      //
949      equalTowZN1LG[gi] = ZN1ADCCorrLG[gi]*equalCoeffZN1[gi];
950      equalTowZP1LG[gi] = ZP1ADCCorrLG[gi]*equalCoeffZP1[gi];
951      equalTowZN2LG[gi] = ZN2ADCCorrLG[gi]*equalCoeffZN2[gi];
952      equalTowZP2LG[gi] = ZP2ADCCorrLG[gi]*equalCoeffZP2[gi];
953   }
954   
955   // Energy calibration of detector responses
956   Float_t calibTowZN1HG[5], calibTowZN2HG[5], calibTowZP1HG[5], calibTowZP2HG[5];
957   Float_t calibSumZN1HG=0., calibSumZN2HG=0., calibSumZP1HG=0., calibSumZP2HG=0.;
958   Float_t calibTowZN1LG[5], calibTowZN2LG[5], calibTowZP1LG[5], calibTowZP2LG[5];
959   Float_t calibSumZN1LG=0., calibSumZN2LG=0., calibSumZ12LG=0., calibSumZP2LG=0.;
960   for(Int_t gi=0; gi<5; gi++){
961      calibTowZN1HG[gi] = equalTowZN1HG[gi]*calibEne[0];
962      calibTowZP1HG[gi] = equalTowZP1HG[gi]*calibEne[1];
963      calibTowZN2HG[gi] = equalTowZN2HG[gi]*calibEne[2];
964      calibTowZP2HG[gi] = equalTowZP2HG[gi]*calibEne[3];
965      calibSumZN1HG += calibTowZN1HG[gi];
966      calibSumZP1HG += calibTowZP1HG[gi];
967      calibSumZN2HG += calibTowZN2HG[gi];
968      calibSumZP2HG += calibTowZP2HG[gi];
969      //
970      calibTowZN1LG[gi] = equalTowZN1LG[gi]*calibEne[0];
971      calibTowZP1LG[gi] = equalTowZP1LG[gi]*calibEne[1];
972      calibTowZN2LG[gi] = equalTowZN2LG[gi]*calibEne[2];
973      calibTowZP2LG[gi] = equalTowZP2LG[gi]*calibEne[3];
974      calibSumZN1LG += calibTowZN1LG[gi];
975      calibSumZ12LG += calibTowZP1LG[gi];
976      calibSumZN2LG += calibTowZN2LG[gi];
977      calibSumZP2LG += calibTowZP2LG[gi];
978   }
979   
980   //  ---      Number of detected spectator nucleons
981   //  *** N.B. -> It works only in Pb-Pb
982   Int_t nDetSpecNLeft, nDetSpecPLeft, nDetSpecNRight, nDetSpecPRight;
983   nDetSpecNLeft = (Int_t) (calibSumZN1HG/2.760);
984   nDetSpecPLeft = (Int_t) (calibSumZP1HG/2.760);
985   nDetSpecNRight = (Int_t) (calibSumZN2HG/2.760);
986   nDetSpecPRight = (Int_t) (calibSumZP2HG/2.760);
987   /*printf("\n\t AliZDCReconstructor -> nDetSpecNLeft %d, nDetSpecPLeft %d,"
988     " nDetSpecNRight %d, nDetSpecPRight %d\n",nDetSpecNLeft, nDetSpecPLeft, 
989     nDetSpecNRight, nDetSpecPRight);*/
990
991   //  ---      Number of generated spectator nucleons (from HIJING parameterization)
992   Int_t nGenSpecNLeft=0, nGenSpecPLeft=0, nGenSpecLeft=0;
993   Int_t nGenSpecNRight=0, nGenSpecPRight=0, nGenSpecRight=0;
994   Double_t impPar=0.;
995   //
996   // *** RECONSTRUCTION FROM SIMULATED DATA
997   // Cut value for Ezem (GeV)
998   // ### Results from production  -> 0<b<18 fm (Apr 2002)
999   /*Float_t eZEMCut = 420.;
1000   Float_t deltaEZEMSup = 690.; 
1001   Float_t deltaEZEMInf = 270.; 
1002   if(zemenergy > (eZEMCut+deltaEZEMSup)){
1003     nGenSpecNLeft  = (Int_t) (fZNCen->Eval(ZN1CalibSum));
1004     nGenSpecPLeft  = (Int_t) (fZPCen->Eval(ZP1CalibSum));
1005     nGenSpecLeft   = (Int_t) (fZDCCen->Eval(ZN1CalibSum+ZP1CalibSum));
1006     nGenSpecNRight = (Int_t) (fZNCen->Eval(ZN2CalibSum));
1007     nGenSpecPRight = (Int_t) (fZNCen->Eval(ZP2CalibSum));
1008     nGenSpecRight  = (Int_t) (fZNCen->Eval(ZN2CalibSum+ZP2CalibSum));
1009     impPar  = fbCen->Eval(ZN1CalibSum+ZP1CalibSum);
1010   }
1011   else if(zemenergy < (eZEMCut-deltaEZEMInf)){
1012     nGenSpecNLeft = (Int_t) (fZNPer->Eval(ZN1CalibSum)); 
1013     nGenSpecPLeft = (Int_t) (fZPPer->Eval(ZP1CalibSum));
1014     nGenSpecLeft  = (Int_t) (fZDCPer->Eval(ZN1CalibSum+ZP1CalibSum));
1015     impPar   = fbPer->Eval(ZN1CalibSum+ZP1CalibSum);
1016   }
1017   else if(zemenergy >= (eZEMCut-deltaEZEMInf) && zemenergy <= (eZEMCut+deltaEZEMSup)){
1018     nGenSpecNLeft = (Int_t) (fZEMn->Eval(zemenergy));
1019     nGenSpecPLeft = (Int_t) (fZEMp->Eval(zemenergy));
1020     nGenSpecLeft  = (Int_t)(fZEMsp->Eval(zemenergy));
1021     impPar   =  fZEMb->Eval(zemenergy);
1022   }
1023   // ### Results from production  -> 0<b<18 fm (Apr 2002)
1024   if(ZN1CalibSum>162.)  nGenSpecNLeft = (Int_t) (fZEMn->Eval(zemenergy));
1025   if(ZP1CalibSum>59.75)  nGenSpecPLeft = (Int_t) (fZEMp->Eval(zemenergy));
1026   if(ZN1CalibSum+ZP1CalibSum>221.5) nGenSpecLeft  = (Int_t)(fZEMsp->Eval(zemenergy));
1027   if(ZN1CalibSum+ZP1CalibSum>220.)  impPar    =  fZEMb->Eval(zemenergy);
1028   */
1029   //
1030   //
1031   // *** RECONSTRUCTION FROM REAL DATA
1032   //
1033   if(corrADCZEMHG > supValueZEM){
1034     nGenSpecNLeft  = (Int_t) (fZNCen->Eval(calibSumZN1HG));
1035     nGenSpecPLeft  = (Int_t) (fZPCen->Eval(calibSumZP1HG));
1036     nGenSpecLeft   = (Int_t) (fZDCCen->Eval(calibSumZN1HG+calibSumZP1HG));
1037     nGenSpecNRight = (Int_t) (fZNCen->Eval(calibSumZN2HG));
1038     nGenSpecPRight = (Int_t) (fZNCen->Eval(calibSumZP2HG));
1039     nGenSpecRight  = (Int_t) (fZNCen->Eval(calibSumZN2HG+calibSumZP2HG));
1040     impPar  = fbCen->Eval(calibSumZN1HG+calibSumZP1HG);
1041   }
1042   else if(corrADCZEMHG < infValueZEM){
1043     nGenSpecNLeft = (Int_t) (fZNPer->Eval(calibSumZN1HG)); 
1044     nGenSpecPLeft = (Int_t) (fZPPer->Eval(calibSumZP1HG));
1045     nGenSpecLeft  = (Int_t) (fZDCPer->Eval(calibSumZN1HG+calibSumZP1HG));
1046     impPar   = fbPer->Eval(calibSumZN1HG+calibSumZP1HG);
1047   }
1048   else if(corrADCZEMHG >= infValueZEM && corrADCZEMHG <= supValueZEM){
1049     nGenSpecNLeft = (Int_t) (fZEMn->Eval(corrADCZEMHG));
1050     nGenSpecPLeft = (Int_t) (fZEMp->Eval(corrADCZEMHG));
1051     nGenSpecLeft  = (Int_t)(fZEMsp->Eval(corrADCZEMHG));
1052     impPar   =  fZEMb->Eval(corrADCZEMHG);
1053   }
1054   // 
1055   if(calibSumZN1HG/maxValEZN1>1.)  nGenSpecNLeft = (Int_t) (fZEMn->Eval(corrADCZEMHG));
1056   if(calibSumZP1HG/maxValEZP1>1.)  nGenSpecPLeft = (Int_t) (fZEMp->Eval(corrADCZEMHG));
1057   if((calibSumZN1HG+calibSumZP1HG/maxValEZDC1)>1.){
1058      nGenSpecLeft = (Int_t)(fZEMsp->Eval(corrADCZEMHG));
1059      impPar = fZEMb->Eval(corrADCZEMHG);
1060   }
1061   if(calibSumZN2HG/maxValEZN2>1.)  nGenSpecNRight = (Int_t) (fZEMn->Eval(corrADCZEMHG));
1062   if(calibSumZP2HG/maxValEZP2>1.)  nGenSpecPRight = (Int_t) (fZEMp->Eval(corrADCZEMHG));
1063   if((calibSumZN2HG+calibSumZP2HG/maxValEZDC2)>1.) nGenSpecRight = (Int_t)(fZEMsp->Eval(corrADCZEMHG));
1064   //
1065   if(nGenSpecNLeft>125)    nGenSpecNLeft=125;
1066   else if(nGenSpecNLeft<0) nGenSpecNLeft=0;
1067   if(nGenSpecPLeft>82)     nGenSpecPLeft=82;
1068   else if(nGenSpecPLeft<0) nGenSpecPLeft=0;
1069   if(nGenSpecLeft>207)     nGenSpecLeft=207;
1070   else if(nGenSpecLeft<0)  nGenSpecLeft=0;
1071   
1072   //  ---      Number of generated participants (from HIJING parameterization)
1073   Int_t nPart, nPartTotLeft, nPartTotRight;
1074   nPart = 207-nGenSpecNLeft-nGenSpecPLeft;
1075   nPartTotLeft = 207-nGenSpecLeft;
1076   nPartTotRight = 207-nGenSpecRight;
1077   if(nPart<0) nPart=0;
1078   if(nPartTotLeft<0) nPartTotLeft=0;
1079   if(nPartTotRight<0) nPartTotRight=0;
1080   //
1081   // *** DEBUG ***
1082   printf("\n\t AliZDCReconstructor -> calibSumZN1HG %1.0f, calibSumZP1HG %1.0f,"
1083       "  calibSumZN2HG %1.0f, calibSumZP2HG %1.0f, corrADCZEMHG %1.0f\n", 
1084       calibSumZN1HG,calibSumZP1HG,calibSumZN2HG,calibSumZP2HG,corrADCZEMHG);
1085   printf("\t AliZDCReconstructor -> nGenSpecNLeft %d, nGenSpecPLeft %d, nGenSpecLeft %d\n"
1086       "\t\t nGenSpecNRight %d, nGenSpecPRight %d, nGenSpecRight %d\n", 
1087       nGenSpecNLeft, nGenSpecPLeft, nGenSpecLeft, 
1088       nGenSpecNRight, nGenSpecPRight, nGenSpecRight);
1089   printf("\t AliZDCReconstructor ->  NpartL %d,  NpartR %d,  b %1.2f fm\n\n",nPartTotLeft, nPartTotRight, impPar);
1090
1091   // create the output tree
1092   AliZDCReco reco(calibSumZN1HG, calibSumZP1HG, calibSumZN2HG, calibSumZP2HG, 
1093                   calibTowZN1LG, calibTowZN2LG, calibTowZP1LG, calibTowZP2LG, 
1094                   corrADCZEMHG, 
1095                   nDetSpecNLeft, nDetSpecPLeft, nDetSpecNRight, nDetSpecPRight, 
1096                   nGenSpecNLeft, nGenSpecPLeft, nGenSpecLeft, nGenSpecNRight, 
1097                   nGenSpecPRight, nGenSpecRight,
1098                   nPartTotLeft, nPartTotRight, impPar);
1099                   
1100   AliZDCReco* preco = &reco;
1101   const Int_t kBufferSize = 4000;
1102   clustersTree->Branch("ZDC", "AliZDCReco", &preco, kBufferSize);
1103
1104   // write the output tree
1105   clustersTree->Fill();
1106 }
1107
1108 //_____________________________________________________________________________
1109 void AliZDCReconstructor::FillZDCintoESD(TTree *clustersTree, AliESDEvent* esd) const
1110 {
1111   // fill energies and number of participants to the ESD
1112
1113   AliZDCReco reco;
1114   AliZDCReco* preco = &reco;
1115   clustersTree->SetBranchAddress("ZDC", &preco);
1116
1117   clustersTree->GetEntry(0);
1118   /*Double_t tZN1Ene[4], tZN2Ene[4];
1119   for(Int_t i=0; i<4; i++){
1120      tZN1Ene[i] = reco.GetZN1EnTow(i);
1121      tZN2Ene[i] = reco.GetZN2EnTow(i);
1122   }
1123   esd->SetZDC(tZN1Ene, tZN2Ene, reco.GetZN1Energy(), reco.GetZP1Energy(), reco.GetZEMsignal(),
1124               reco.GetZN2Energy(), reco.GetZP2Energy(), 
1125               reco.GetNPartLeft());
1126   */
1127   esd->SetZDC(reco.GetZN1Energy(), reco.GetZP1Energy(), reco.GetZEMsignal(),
1128               reco.GetZN2Energy(), reco.GetZP2Energy(), 
1129               reco.GetNPartLeft());
1130   
1131   /*Double_t tZN1Ene[4], tZN2Ene[4];
1132   for(Int_t i=0; i<4; i++){
1133      tZN1Ene[i] = reco.GetZN1EnTow(i);
1134      tZN2Ene[i] = reco.GetZN2EnTow(i);
1135   }*/
1136   //esd->SetZN1TowerEnergy(tZN1Ene);
1137   //esd->SetZN2TowerEnergy(tZN2Ene);
1138   
1139 }
1140
1141 //_____________________________________________________________________________
1142 AliCDBStorage* AliZDCReconstructor::SetStorage(const char *uri) 
1143 {
1144   // Setting the storage
1145
1146   Bool_t deleteManager = kFALSE;
1147   
1148   AliCDBManager *manager = AliCDBManager::Instance();
1149   AliCDBStorage *defstorage = manager->GetDefaultStorage();
1150   
1151   if(!defstorage || !(defstorage->Contains("ZDC"))){ 
1152      AliWarning("No default storage set or default storage doesn't contain ZDC!");
1153      manager->SetDefaultStorage(uri);
1154      deleteManager = kTRUE;
1155   }
1156  
1157   AliCDBStorage *storage = manager->GetDefaultStorage();
1158
1159   if(deleteManager){
1160     AliCDBManager::Instance()->UnsetDefaultStorage();
1161     defstorage = 0;   // the storage is killed by AliCDBManager::Instance()->Destroy()
1162   }
1163
1164   return storage; 
1165 }
1166
1167 //_____________________________________________________________________________
1168 AliZDCCalibData* AliZDCReconstructor::GetCalibData() const
1169 {
1170
1171   // Getting calibration object for ZDC set
1172
1173   AliCDBEntry  *entry = AliCDBManager::Instance()->Get("ZDC/Calib/Data");
1174   if(!entry) AliFatal("No calibration data loaded!");  
1175
1176   AliZDCCalibData *calibdata = dynamic_cast<AliZDCCalibData*>  (entry->GetObject());
1177   if(!calibdata)  AliFatal("Wrong calibration object in calibration  file!");
1178
1179   return calibdata;
1180 }