Correct energy calibration in case of p-p data
[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 //                  Author: Chiara.Oppedisano@to.infn.it                     //
22 //                                                                           //
23 // NOTATIONS ADOPTED TO IDENTIFY DETECTORS (used in different ages!):        //
24 //   (ZN1,ZP1) or (ZNC, ZPC) or RIGHT refers to side C (RB26)                //
25 //   (ZN2,ZP2) or (ZNA, ZPA) or LEFT refers to side A (RB24)                 //
26 //                                                                           //
27 ///////////////////////////////////////////////////////////////////////////////
28
29
30 #include <TH2F.h>
31 #include <TH1D.h>
32 #include <TAxis.h>
33 #include <TMap.h>
34
35 #include "AliRawReader.h"
36 #include "AliGRPObject.h"
37 #include "AliESDEvent.h"
38 #include "AliESDZDC.h"
39 #include "AliZDCDigit.h"
40 #include "AliZDCRawStream.h"
41 #include "AliZDCReco.h"
42 #include "AliZDCReconstructor.h"
43 #include "AliZDCPedestals.h"
44 #include "AliZDCEnCalib.h"
45 #include "AliZDCTowerCalib.h"
46 #include "AliZDCRecoParam.h"
47 #include "AliZDCRecoParampp.h"
48 #include "AliZDCRecoParamPbPb.h"
49 #include "AliRunInfo.h"
50
51
52 ClassImp(AliZDCReconstructor)
53 AliZDCRecoParam *AliZDCReconstructor::fRecoParam=0;  //reconstruction parameters
54
55 //_____________________________________________________________________________
56 AliZDCReconstructor:: AliZDCReconstructor() :
57   fPedData(GetPedestalData()),
58   fEnCalibData(GetEnergyCalibData()),
59   fTowCalibData(GetTowerCalibData()),
60   fRecoMode(0),
61   fBeamEnergy(0.),
62   fNRun(0),
63   fIsCalibrationMB(kFALSE),
64   fPedSubMode(0),
65   fSignalThreshold(7)
66 {
67   // **** Default constructor
68
69 }
70
71
72 //_____________________________________________________________________________
73 AliZDCReconstructor::~AliZDCReconstructor()
74 {
75 // destructor
76 //   if(fRecoParam)    delete fRecoParam;
77    if(fPedData)      delete fPedData;    
78    if(fEnCalibData)  delete fEnCalibData;
79    if(fTowCalibData) delete fTowCalibData;
80 }
81
82 //____________________________________________________________________________
83 void AliZDCReconstructor::Init()
84 {
85   // Setting reconstruction mode
86   // Getting beam type and beam energy from GRP calibration object
87   
88   TString runType = GetRunInfo()->GetRunType();
89
90   if((runType.CompareTo("CALIBRATION_MB")) == 0){
91     fIsCalibrationMB = kTRUE;
92   }
93     
94   TString beamType = GetRunInfo()->GetBeamType();
95   // This is a temporary solution to allow reconstruction in tests without beam
96   if(((beamType.CompareTo("UNKNOWN"))==0) && ((runType.CompareTo("PHYSICS")) == 0)){
97     fRecoMode=1;
98   }
99   else if(beamType==AliGRPObject::GetInvalidString()){
100     AliWarning("GRP/GRP/Data entry:  missing value for the beam type !");
101     AliError("\t ZDC does not reconstruct event 4 UNKNOWN beam type\n");
102     return;
103   }
104
105   if(((beamType.CompareTo("pp"))==0) || ((beamType.CompareTo("p-p"))==0)
106      ||((beamType.CompareTo("PP"))==0) || ((beamType.CompareTo("P-P"))==0)) fRecoMode=1;
107   else if((beamType.CompareTo("A-A")) == 0){
108     fRecoMode=2;
109     if(fIsCalibrationMB == kTRUE){ 
110       fRecoParam = new AliZDCRecoParamPbPb();
111       //
112       TH2F* hZDCvsZEM = new TH2F("hZDCvsZEM","hZDCvsZEM",100,0.,10.,100,0.,1000.);
113       hZDCvsZEM->SetXTitle("E_{ZEM} (TeV)"); hZDCvsZEM->SetYTitle("E_{ZDC} (TeV)");
114       fRecoParam->SetZDCvsZEM(hZDCvsZEM);
115       //
116       TH2F* hZDCCvsZEM = new TH2F("hZDCCvsZEM","hZDCCvsZEM",100,0.,10.,100,0.,500.);
117       hZDCCvsZEM->SetXTitle("E_{ZEM} (TeV)"); hZDCCvsZEM->SetYTitle("E_{ZDCC} (TeV)");
118       fRecoParam->SetZDCCvsZEM(hZDCCvsZEM);
119       //
120       TH2F* hZDCAvsZEM = new TH2F("hZDCAvsZEM","hZDCAvsZEM",100,0.,10.,100,0.,500.);
121       hZDCAvsZEM->SetXTitle("E_{ZEM} (TeV)"); hZDCAvsZEM->SetYTitle("E_{ZDCA} (TeV)"); 
122       fRecoParam->SetZDCAvsZEM(hZDCAvsZEM);
123       //
124       AliInfo("\n ***** CALIBRATION_MB data -> building AliZDCRecoParamPbPb object *****");
125     }
126   }
127     
128   fBeamEnergy = GetRunInfo()->GetBeamEnergy();
129   if(fBeamEnergy==AliGRPObject::GetInvalidFloat()){
130     AliWarning("GRP/GRP/Data entry:  missing value for the beam energy ! Using 0.");
131     fBeamEnergy = 0.;
132   }
133     
134   if(fIsCalibrationMB==kFALSE)  
135     printf("\n\n ***** ZDC reconstruction initialized for %s @ %1.0f GeV *****\n\n",beamType.Data(), fBeamEnergy);
136   
137 }
138
139 //_____________________________________________________________________________
140 void AliZDCReconstructor::Reconstruct(TTree* digitsTree, TTree* clustersTree) const
141 {
142   // *** Local ZDC reconstruction for digits
143   // Works on the current event
144     
145   // Retrieving calibration data  
146   // Parameters for mean value pedestal subtraction
147   int const kNch = 24;
148   Float_t meanPed[2*kNch];    
149   for(Int_t jj=0; jj<2*kNch; jj++) meanPed[jj] = fPedData->GetMeanPed(jj);
150   // Parameters pedestal subtraction through correlation with out-of-time signals
151   Float_t corrCoeff0[2*kNch], corrCoeff1[2*kNch];
152   for(Int_t jj=0; jj<2*kNch; jj++){
153      corrCoeff0[jj] =  fPedData->GetPedCorrCoeff0(jj);
154      corrCoeff1[jj] =  fPedData->GetPedCorrCoeff1(jj);
155   }
156
157   // get digits
158   AliZDCDigit digit;
159   AliZDCDigit* pdigit = &digit;
160   digitsTree->SetBranchAddress("ZDC", &pdigit);
161   //printf("\n\t # of digits in tree: %d\n",(Int_t) digitsTree->GetEntries());
162
163   // loop over digits
164   Float_t tZN1Corr[10], tZP1Corr[10], tZN2Corr[10], tZP2Corr[10]; 
165   Float_t dZEM1Corr[2], dZEM2Corr[2], sPMRef1[2], sPMRef2[2]; 
166   for(Int_t i=0; i<10; i++){
167      tZN1Corr[i] = tZP1Corr[i] = tZN2Corr[i] = tZP2Corr[i] = 0.;
168      if(i<2) dZEM1Corr[i] = dZEM2Corr[i] = sPMRef1[i] = sPMRef2[i] = 0.;
169   }  
170   
171   Int_t digNentries = digitsTree->GetEntries();
172   Float_t ootDigi[kNch];
173   // -- Reading out-of-time signals (last kNch entries) for current event
174   if(fPedSubMode==1){
175     for(Int_t iDigit=kNch; iDigit<digNentries; iDigit++){
176        ootDigi[iDigit] = digitsTree->GetEntry(iDigit);
177     }
178   }
179   
180   for(Int_t iDigit=0; iDigit<(digNentries/2); iDigit++) {
181    digitsTree->GetEntry(iDigit);
182    if (!pdigit) continue;
183    //  
184    Int_t det = digit.GetSector(0);
185    Int_t quad = digit.GetSector(1);
186    Int_t pedindex = -1;
187    Float_t ped2SubHg=0., ped2SubLg=0.;
188    if(quad!=5){
189      if(det==1)      pedindex = quad;
190      else if(det==2) pedindex = quad+5;
191      else if(det==3) pedindex = quad+9;
192      else if(det==4) pedindex = quad+12;
193      else if(det==5) pedindex = quad+17;
194    }
195    else pedindex = (det-1)/3+22;
196    //
197    if(fPedSubMode==0){
198      ped2SubHg = meanPed[pedindex];
199      ped2SubLg = meanPed[pedindex+kNch];
200    }
201    else if(fPedSubMode==1){
202      ped2SubHg = corrCoeff1[pedindex]*ootDigi[pedindex]+corrCoeff0[pedindex];
203      ped2SubLg = corrCoeff1[pedindex+kNch]*ootDigi[pedindex+kNch]+corrCoeff0[pedindex+kNch];
204    }
205       
206    if(quad != 5){ // ZDC (not reference PTMs!)
207     if(det == 1){ // *** ZNC
208        tZN1Corr[quad] = (Float_t) (digit.GetADCValue(0)-ped2SubHg);
209        tZN1Corr[quad+5] = (Float_t) (digit.GetADCValue(1)-ped2SubLg);
210        if(tZN1Corr[quad]<0.) tZN1Corr[quad] = 0.;
211        if(tZN1Corr[quad+5]<0.) tZN1Corr[quad+5] = 0.;
212     }
213     else if(det == 2){ // *** ZP1
214        tZP1Corr[quad] = (Float_t) (digit.GetADCValue(0)-ped2SubHg);
215        tZP1Corr[quad+5] = (Float_t) (digit.GetADCValue(1)-ped2SubLg);
216        if(tZP1Corr[quad]<0.) tZP1Corr[quad] = 0.;
217        if(tZP1Corr[quad+5]<0.) tZP1Corr[quad+5] = 0.;
218     }
219     else if(det == 3){
220        if(quad == 1){       // *** ZEM1  
221          dZEM1Corr[0] += (Float_t) (digit.GetADCValue(0)-ped2SubHg); 
222          dZEM1Corr[1] += (Float_t) (digit.GetADCValue(1)-ped2SubLg); 
223          if(dZEM1Corr[0]<0.) dZEM1Corr[0] = 0.;
224          if(dZEM1Corr[1]<0.) dZEM1Corr[1] = 0.;
225        }
226        else if(quad == 2){  // *** ZEM2
227          dZEM2Corr[0] += (Float_t) (digit.GetADCValue(0)-ped2SubHg); 
228          dZEM2Corr[1] += (Float_t) (digit.GetADCValue(1)-ped2SubLg); 
229          if(dZEM2Corr[0]<0.) dZEM2Corr[0] = 0.;
230          if(dZEM2Corr[1]<0.) dZEM2Corr[1] = 0.;
231        }
232     }
233     else if(det == 4){  // *** ZN2
234        tZN2Corr[quad] = (Float_t) (digit.GetADCValue(0)-ped2SubHg);
235        tZN2Corr[quad+5] = (Float_t) (digit.GetADCValue(1)-ped2SubLg);
236        if(tZN2Corr[quad]<0.) tZN2Corr[quad] = 0.;
237        if(tZN2Corr[quad+5]<0.) tZN2Corr[quad+5] = 0.;
238    }
239     else if(det == 5){  // *** ZP2 
240        tZP2Corr[quad] = (Float_t) (digit.GetADCValue(0)-ped2SubHg);
241        tZP2Corr[quad+5] = (Float_t) (digit.GetADCValue(1)-ped2SubLg);
242        if(tZP2Corr[quad]<0.) tZP2Corr[quad] = 0.;
243        if(tZP2Corr[quad+5]<0.) tZP2Corr[quad+5] = 0.;
244     }
245    }
246    else{ // Reference PMs
247      if(det == 1){
248        sPMRef1[0] = (Float_t) (digit.GetADCValue(0)-ped2SubHg);
249        sPMRef1[1] = (Float_t) (digit.GetADCValue(1)-ped2SubLg);
250        // Ch. debug
251        if(sPMRef1[0]<0.) sPMRef1[0] = 0.;
252        if(sPMRef2[1]<0.) sPMRef1[1] = 0.;
253      }
254      else if(det == 4){
255        sPMRef2[0] = (Float_t) (digit.GetADCValue(0)-ped2SubHg);
256        sPMRef2[1] = (Float_t) (digit.GetADCValue(1)-ped2SubLg);
257        // Ch. debug
258        if(sPMRef2[0]<0.) sPMRef2[0] = 0.;
259        if(sPMRef2[1]<0.) sPMRef2[1] = 0.;
260      }
261    }
262
263    // Ch. debug
264    /*printf("AliZDCReconstructor: digit #%d det %d quad %d pedHG %1.0f pedLG %1.0f\n",
265          iDigit, det, quad, ped2SubHg, ped2SubLg);
266    printf(" -> pedindex %d\n", pedindex);
267    printf("   HGChain -> RawDig %d DigCorr %1.2f", 
268         digit.GetADCValue(0), digit.GetADCValue(0)-ped2SubHg); 
269    printf("   LGChain -> RawDig %d DigCorr %1.2f\n", 
270         digit.GetADCValue(1), digit.GetADCValue(1)-ped2SubLg);*/ 
271    
272   }//digits loop
273  
274   
275   // If CALIBRATION_MB run build the RecoParam object 
276   if(fIsCalibrationMB){
277     Float_t ZDCC=0., ZDCA=0., ZEM=0;
278     ZEM += dZEM1Corr[0] + dZEM2Corr[0];
279     for(Int_t jkl=0; jkl<5; jkl++){
280        ZDCC += tZN1Corr[jkl] + tZP1Corr[jkl];
281        ZDCA += tZN2Corr[jkl] + tZP2Corr[jkl];
282     }
283     // Using energies in TeV in fRecoParam object
284     BuildRecoParam(ZDCC/1000., ZDCA/1000., ZEM/1000.);
285   }
286   Bool_t recFlag1 = kFALSE, recFlag2 = kFALSE, recFlag3 = kFALSE;
287   // reconstruct the event
288   if(fRecoMode==1)
289     ReconstructEventpp(clustersTree, tZN1Corr, tZP1Corr, tZN2Corr, tZP2Corr, 
290       dZEM1Corr, dZEM2Corr, sPMRef1, sPMRef2, recFlag1, recFlag2, recFlag3);
291   else if(fRecoMode==2)
292     ReconstructEventPbPb(clustersTree, tZN1Corr, tZP1Corr, tZN2Corr, tZP2Corr, 
293       dZEM1Corr, dZEM2Corr, sPMRef1, sPMRef2, recFlag1, recFlag2, recFlag3);
294 }
295
296 //_____________________________________________________________________________
297 void AliZDCReconstructor::Reconstruct(AliRawReader* rawReader, TTree* clustersTree) const
298 {
299   // *** ZDC raw data reconstruction
300   // Works on the current event
301   
302   // Retrieving calibration data  
303   // Parameters for pedestal subtraction
304   int const kNch = 24;
305   Float_t meanPed[2*kNch];    
306   for(Int_t jj=0; jj<2*kNch; jj++) meanPed[jj] = fPedData->GetMeanPed(jj);
307   // Parameters pedestal subtraction through correlation with out-of-time signals
308   Float_t corrCoeff0[2*kNch], corrCoeff1[2*kNch];
309   for(Int_t jj=0; jj<2*kNch; jj++){
310      corrCoeff0[jj] =  fPedData->GetPedCorrCoeff0(jj);
311      corrCoeff1[jj] =  fPedData->GetPedCorrCoeff1(jj);
312   }
313
314   Int_t adcZN1[5], adcZN1oot[5], adcZN1lg[5], adcZN1ootlg[5];
315   Int_t adcZP1[5], adcZP1oot[5], adcZP1lg[5], adcZP1ootlg[5];
316   Int_t adcZN2[5], adcZN2oot[5], adcZN2lg[5], adcZN2ootlg[5];
317   Int_t adcZP2[5], adcZP2oot[5], adcZP2lg[5], adcZP2ootlg[5];
318   Int_t adcZEM[2], adcZEMoot[2], adcZEMlg[2], adcZEMootlg[2];
319   Int_t pmRef[2], pmRefoot[2], pmReflg[2], pmRefootlg[2];
320   for(Int_t ich=0; ich<5; ich++){
321     adcZN1[ich] = adcZN1oot[ich] = adcZN1lg[ich] = adcZN1ootlg[ich] = 0;
322     adcZP1[ich] = adcZP1oot[ich] = adcZP1lg[ich] = adcZP1ootlg[ich] = 0;
323     adcZN2[ich] = adcZN2oot[ich] = adcZN2lg[ich] = adcZN2ootlg[ich] = 0;
324     adcZP2[ich] = adcZP2oot[ich] = adcZP2lg[ich] = adcZP2ootlg[ich] = 0;
325     if(ich<2){
326       adcZEM[ich] = adcZEMoot[ich] = adcZEMlg[ich] = adcZEMootlg[ich] = 0;
327       pmRef[ich] = pmRefoot[ich] = pmReflg[ich] = pmRefootlg[ich] = 0;
328     }
329   }
330   
331   // loop over raw data
332   Float_t tZN1Corr[10], tZP1Corr[10], tZN2Corr[10], tZP2Corr[10]; 
333   Float_t dZEM1Corr[2], dZEM2Corr[2], sPMRef1[2], sPMRef2[2]; 
334   for(Int_t i=0; i<10; i++){
335      tZN1Corr[i] = tZP1Corr[i] = tZN2Corr[i] = tZP2Corr[i] = 0.;
336      if(i<2) dZEM1Corr[i] = dZEM2Corr[i] = sPMRef1[i] = sPMRef2[i] = 0.;
337   }  
338
339   //fNRun = (Int_t) rawReader->GetRunNumber();
340   Bool_t chOff=kFALSE, isUndflw=kFALSE, isOvflw=kFALSE, isADCEvGood=kFALSE;
341   Int_t kFirstADCGeo=0, kLastADCGeo=3;
342   //Int_t kScalerGeo=8, kPUGeo=29, kTrigScales=30, kTrigHistory=31;
343   //
344   rawReader->Reset();
345   AliZDCRawStream rawData(rawReader);
346   while(rawData.Next()){
347    //if(rawData.IsCalibration() == kFALSE){ // Reading ADCs
348    if((rawData.GetADCModule()>=kFirstADCGeo) && (rawData.GetADCModule()<=kLastADCGeo)){    
349     //printf(" **** Reading ADC raw data from module %d **** \n",rawData.GetADCModule());
350     //
351     if((rawData.IsADCDataWord()) && (rawData.GetNChannelsOn()<48))  chOff = kTRUE;
352     if((rawData.IsADCDataWord()) && (rawData.IsOverflow() == kTRUE)) isUndflw = kTRUE;
353     if((rawData.IsADCDataWord()) && (rawData.IsUnderflow() == kTRUE)) isOvflw = kTRUE;
354     if((rawData.IsADCDataWord()) && (rawData.IsADCEventGood() == kTRUE)) isADCEvGood = kTRUE;
355     
356     if((rawData.IsADCDataWord()) && (isUndflw==kFALSE) 
357         && (isOvflw==kFALSE) && (isADCEvGood==kTRUE)){
358      
359       Int_t adcMod = rawData.GetADCModule();
360       Int_t det = rawData.GetSector(0);
361       Int_t quad = rawData.GetSector(1);
362       Int_t gain = rawData.GetADCGain();
363       Int_t pedindex=0;
364       //
365       // Mean pedestal value subtraction -------------------------------------------------------
366       if(fPedSubMode == 0){
367        // Not interested in o.o.t. signals (ADC modules 2, 3)
368        if(adcMod == 2 || adcMod == 3) continue;
369        //
370        if(quad != 5){ // ZDCs (not reference PTMs)
371         if(det == 1){    
372           pedindex = quad;
373           if(gain == 0) tZN1Corr[quad]  += (Float_t) (rawData.GetADCValue()-meanPed[pedindex]); 
374           else tZN1Corr[quad+5]  += (Float_t) (rawData.GetADCValue()-meanPed[pedindex+kNch]); 
375         }
376         else if(det == 2){ 
377           pedindex = quad+5;
378           if(gain == 0) tZP1Corr[quad]  += (Float_t) (rawData.GetADCValue()-meanPed[pedindex]); 
379           else tZP1Corr[quad+5]  += (Float_t) (rawData.GetADCValue()-meanPed[pedindex+kNch]); 
380         }
381         else if(det == 3){ 
382           pedindex = quad+9;
383           if(quad==1){     
384             if(gain == 0) dZEM1Corr[0] += (Float_t) (rawData.GetADCValue()-meanPed[pedindex]); 
385             else dZEM1Corr[1] += (Float_t) (rawData.GetADCValue()-meanPed[pedindex+kNch]); 
386           }
387           else if(quad==2){ 
388             if(gain == 0) dZEM2Corr[0] += (Float_t) (rawData.GetADCValue()-meanPed[pedindex]); 
389             else dZEM2Corr[1] += (Float_t) (rawData.GetADCValue()-meanPed[pedindex+kNch]); 
390           }
391         }
392         else if(det == 4){       
393           pedindex = quad+12;
394           if(gain == 0) tZN2Corr[quad]  += (Float_t) (rawData.GetADCValue()-meanPed[pedindex]); 
395           else tZN2Corr[quad+5]  += (Float_t) (rawData.GetADCValue()-meanPed[pedindex+kNch]); 
396         }
397         else if(det == 5){
398           pedindex = quad+17;
399           if(gain == 0) tZP2Corr[quad]  += (Float_t) (rawData.GetADCValue()-meanPed[pedindex]); 
400           else tZP2Corr[quad+5]  += (Float_t) (rawData.GetADCValue()-meanPed[pedindex+kNch]); 
401         }
402        }
403        else{ // reference PM
404          pedindex = (det-1)/3 + 22;
405          if(det == 1){
406            if(gain==0) sPMRef1[0] += (Float_t) (rawData.GetADCValue()-meanPed[pedindex]);
407          else sPMRef1[1] += (Float_t) (rawData.GetADCValue()-meanPed[pedindex]);
408          }
409          else if(det == 4){
410            if(gain==0) sPMRef2[0] += (Float_t) (rawData.GetADCValue()-meanPed[pedindex]);
411            else sPMRef2[1] += (Float_t) (rawData.GetADCValue()-meanPed[pedindex]);
412          }
413        }
414        // Ch. debug
415        /*printf(" AliZDCReconstructor: det %d quad %d res %d -> Pedestal[%d] %1.0f\n", 
416          det,quad,gain, pedindex, meanPed[pedindex]);
417        printf("                      RawADC %d ADCCorr %1.0f\n", 
418          rawData.GetADCValue(), rawData.GetADCValue()-meanPed[pedindex]);*/
419          
420       }// mean pedestal subtraction
421       // Pedestal subtraction from correlation ------------------------------------------------
422       else if(fPedSubMode == 1){
423        // In time signals
424        if(adcMod==0 || adcMod==1){
425          if(quad != 5){ // signals from ZDCs
426            if(det == 1){
427              if(gain==0) adcZN1[quad] = rawData.GetADCValue();
428              else adcZN1lg[quad] = rawData.GetADCValue();
429            }
430            else if(det == 2){
431              if(gain==0) adcZP1[quad] = rawData.GetADCValue();
432              else adcZP1lg[quad] = rawData.GetADCValue();
433            }
434            else if(det == 3){
435              if(gain==0) adcZEM[quad-1] = rawData.GetADCValue();
436              else adcZEMlg[quad-1] = rawData.GetADCValue();
437            }
438            else if(det == 4){
439              if(gain==0) adcZN2[quad] = rawData.GetADCValue();
440              else adcZN2lg[quad] = rawData.GetADCValue();
441            }
442            else if(det == 5){
443              if(gain==0) adcZP2[quad] = rawData.GetADCValue();
444              else adcZP2lg[quad] = rawData.GetADCValue();
445            }
446          }
447          else{ // signals from reference PM
448             if(gain==0) pmRef[quad-1] = rawData.GetADCValue();
449             else pmReflg[quad-1] = rawData.GetADCValue();
450          }
451        }
452        // Out-of-time pedestals
453        else if(adcMod==2 || adcMod==3){
454          if(quad != 5){ // signals from ZDCs
455            if(det == 1){
456              if(gain==0) adcZN1oot[quad] = rawData.GetADCValue();
457              else adcZN1ootlg[quad] = rawData.GetADCValue();
458            }
459            else if(det == 2){
460              if(gain==0) adcZP1oot[quad] = rawData.GetADCValue();
461              else adcZP1ootlg[quad] = rawData.GetADCValue();
462            }
463            else if(det == 3){
464              if(gain==0) adcZEMoot[quad-1] = rawData.GetADCValue();
465              else adcZEMootlg[quad-1] = rawData.GetADCValue();
466            }
467            else if(det == 4){
468              if(gain==0) adcZN2oot[quad] = rawData.GetADCValue();
469              else adcZN2ootlg[quad] = rawData.GetADCValue();
470            }
471            else if(det == 5){
472              if(gain==0) adcZP2oot[quad] = rawData.GetADCValue();
473              else adcZP2ootlg[quad] = rawData.GetADCValue();
474            }
475          }
476          else{ // signals from reference PM
477             if(gain==0) pmRefoot[quad-1] = rawData.GetADCValue();
478             else pmRefootlg[quad-1] = rawData.GetADCValue();
479          }
480        }
481       } // pedestal subtraction from correlation
482       // Ch. debug
483       /*printf("\t AliZDCReconstructor: det %d quad %d res %d -> Ped[%d] = %1.0f\n", 
484         det,quad,gain, pedindex, meanPed[pedindex]);*/
485     }//IsADCDataWord
486    }// Not raw data from calibration run!
487    else{
488    }
489   }//loop on raw data
490   
491   if(fPedSubMode==1){
492     for(Int_t t=0; t<5; t++){
493        tZN1Corr[t] = adcZN1[t] - (corrCoeff1[t]*adcZN1oot[t]+corrCoeff0[t]);
494        tZN1Corr[t+5] = adcZN1lg[t] - (corrCoeff1[t+kNch]*adcZN1ootlg[t]+corrCoeff0[t+kNch]);
495        //
496        tZP1Corr[t] = adcZP1[t] - (corrCoeff1[t+5]*adcZP1oot[t]+corrCoeff0[t+5]);
497        tZP1Corr[t+5] = adcZP1lg[t] - (corrCoeff1[t+5+kNch]*adcZP1ootlg[t]+corrCoeff0[t+5+kNch]);
498        //
499        tZN2Corr[t] = adcZN2[t] - (corrCoeff1[t+12]*adcZN2oot[t]+corrCoeff0[t+12]);
500        tZN2Corr[t+5] = adcZN2lg[t] - (corrCoeff1[t+12+kNch]*adcZN2ootlg[t]+corrCoeff0[t+12+kNch]);
501        //
502        tZP2Corr[t] = adcZP2[t] - (corrCoeff1[t+17]*adcZP2oot[t]+corrCoeff0[t+17]);
503        tZP2Corr[t+5] = adcZP2lg[t] - (corrCoeff1[t+17+kNch]*adcZP2ootlg[t]+corrCoeff0[t+17+kNch]);
504        // 0---------0 Ch. debug 0---------0
505 /*       printf("\n\n ---------- Debug of pedestal subtraction from correlation ----------\n");
506        printf("\tCorrCoeff0\tCorrCoeff1\n");
507        printf(" ZN1 %d\t%1.0f\t%1.0f\n",t,corrCoeff0[t],corrCoeff1[t]);
508        printf(" ZN1lg %d\t%1.0f\t%1.0f\n",t+kNch,corrCoeff0[t+kNch],corrCoeff1[t+kNch]);
509        printf(" ZP1 %d\t%1.0f\t%1.0f\n",t+5,corrCoeff0[t+5],corrCoeff1[t+5]);
510        printf(" ZP1lg %d\t%1.0f\t%1.0f\n",t+5+kNch,corrCoeff0[t+5+kNch],corrCoeff1[t+5+kNch]);
511        printf(" ZN2 %d\t%1.0f\t%1.0f\n",t+12,corrCoeff0[t+12],corrCoeff1[t+12]);
512        printf(" ZN2lg %d\t%1.0f\t%1.0f\n",t+12+kNch,corrCoeff0[t+12+kNch],corrCoeff1[t+12+kNch]);
513        printf(" ZP2 %d\t%1.0f\t%1.0f\n",t+17,corrCoeff0[t+17],corrCoeff1[t+17]);
514        printf(" ZP2lg %d\t%1.0f\t%1.0f\n",t+17+kNch,corrCoeff0[t+17+kNch],corrCoeff1[t+17+kNch]);
515        
516        printf("ZN1 ->  rawADC %d\tpedestal%1.2f\tcorrADC%1.2f\n",
517                 adcZN1[t],(corrCoeff1[t]*adcZN1oot[t]+corrCoeff0[t]),tZN1Corr[t]);
518        printf(" lg ->  rawADC %d\tpedestal%1.2f\tcorrADC%1.2f\n",
519                 adcZN1lg[t],(corrCoeff1[t+kNch]*adcZN1ootlg[t]+corrCoeff0[t+kNch]),tZN1Corr[t+5]);
520        //
521        printf("ZP1 ->  rawADC %d\tpedestal%1.2f\tcorrADC%1.2f\n",
522                 adcZP1[t],(corrCoeff1[t+5]*adcZP1oot[t]+corrCoeff0[t+5]),tZP1Corr[t]);
523        printf(" lg ->  rawADC %d\tpedestal%1.2f\tcorrADC%1.2f\n",
524                 adcZP1lg[t],(corrCoeff1[t+5+kNch]*adcZP1ootlg[t]+corrCoeff0[t+5+kNch]),tZP1Corr[t+5]);
525        //
526        printf("ZN2 ->  rawADC %d\tpedestal%1.2f\tcorrADC%1.2f\n",
527                 adcZN2[t],(corrCoeff1[t+12]*adcZN2oot[t]+corrCoeff0[t+12]),tZN2Corr[t]);
528        printf(" lg ->  rawADC %d\tpedestal%1.2f\tcorrADC%1.2f\n",
529                 adcZN2lg[t],(corrCoeff1[t+12+kNch]*adcZN2ootlg[t]+corrCoeff0[t+12+kNch]),tZN2Corr[t+5]);
530        //
531        printf("ZP2 ->  rawADC %d\tpedestal%1.2f\tcorrADC%1.2f\n",
532                 adcZP2[t],(corrCoeff1[t+17]*adcZP2oot[t]+corrCoeff0[t+17]),tZP2Corr[t]);
533        printf(" lg ->  rawADC %d\tpedestal%1.2f\tcorrADC%1.2f\n",
534                 adcZP2lg[t],(corrCoeff1[t+17+kNch]*adcZP2ootlg[t]+corrCoeff0[t+17+kNch]),tZP2Corr[t+5]);
535 */
536     }
537     dZEM1Corr[0] = adcZEM[0]   - (corrCoeff1[9]*adcZEMoot[0]+corrCoeff0[9]);
538     dZEM1Corr[1] = adcZEMlg[0] - (corrCoeff1[9+kNch]*adcZEMootlg[0]+corrCoeff0[9+kNch]);
539     dZEM2Corr[0] = adcZEM[1]   - (corrCoeff1[10]*adcZEMoot[1]+corrCoeff0[10]);
540     dZEM2Corr[1] = adcZEMlg[1] - (corrCoeff1[10+kNch]*adcZEMootlg[1]+corrCoeff0[10+kNch]);
541     //
542     sPMRef1[0] = pmRef[0]   - (corrCoeff1[22]*pmRefoot[0]+corrCoeff0[22]);
543     sPMRef1[1] = pmReflg[0] - (corrCoeff1[22+kNch]*pmRefootlg[0]+corrCoeff0[22+kNch]);
544     sPMRef2[0] = pmRef[0]   - (corrCoeff1[23]*pmRefoot[1]+corrCoeff0[23]);
545     sPMRef2[1] = pmReflg[0] - (corrCoeff1[23+kNch]*pmRefootlg[1]+corrCoeff0[23+kNch]);
546   }
547     
548   // If CALIBRATION_MB run build the RecoParam object 
549   if(fIsCalibrationMB){
550     Float_t ZDCC=0., ZDCA=0., ZEM=0;
551     ZEM += dZEM1Corr[0] + dZEM2Corr[0];
552     for(Int_t jkl=0; jkl<5; jkl++){
553        ZDCC += tZN1Corr[jkl] + tZP1Corr[jkl];
554        ZDCA += tZN2Corr[jkl] + tZP2Corr[jkl];
555     }
556     BuildRecoParam(ZDCC/100., ZDCA/100., ZEM/100.);
557   }
558   // reconstruct the event
559   else{
560     if(fRecoMode==1) // p-p data
561       ReconstructEventpp(clustersTree, tZN1Corr, tZP1Corr, tZN2Corr, tZP2Corr, 
562         dZEM1Corr, dZEM2Corr, sPMRef1, sPMRef2, chOff, isUndflw, isOvflw);
563     else if(fRecoMode==2) // Pb-Pb data
564       ReconstructEventPbPb(clustersTree, tZN1Corr, tZP1Corr, tZN2Corr, tZP2Corr, 
565         dZEM1Corr, dZEM2Corr, sPMRef1, sPMRef2, chOff, isUndflw, isOvflw);
566   }
567 }
568
569 //_____________________________________________________________________________
570 void AliZDCReconstructor::ReconstructEventpp(TTree *clustersTree, Float_t* corrADCZN1, 
571         Float_t* corrADCZP1, Float_t* corrADCZN2, Float_t* corrADCZP2,
572         Float_t* corrADCZEM1, Float_t* corrADCZEM2, Float_t* sPMRef1, Float_t* sPMRef2,
573         Bool_t channelsOff, Bool_t chUnderflow, Bool_t chOverflow) const
574 {
575   // ****************** Reconstruct one event ******************
576   
577   // CH. debug
578   /*printf("\n*************************************************\n");
579   printf(" ReconstructEventpp -> values after pedestal subtraction:\n");
580   printf(" ADCZN1 [%1.2f %1.2f %1.2f %1.2f %1.2f]\n",
581         corrADCZN1[0],corrADCZN1[1],corrADCZN1[2],corrADCZN1[3],corrADCZN1[4]);
582   printf(" ADCZP1 [%1.2f %1.2f %1.2f %1.2f %1.2f]\n",
583         corrADCZP1[0],corrADCZP1[1],corrADCZP1[2],corrADCZP1[3],corrADCZP1[4]);
584   printf(" ADCZN2 [%1.2f %1.2f %1.2f %1.2f %1.2f]\n",
585         corrADCZN2[0],corrADCZN2[1],corrADCZN2[2],corrADCZN2[3],corrADCZN2[4]);
586   printf(" ADCZP2 [%1.2f %1.2f %1.2f %1.2f %1.2f]\n",
587         corrADCZP2[0],corrADCZP2[1],corrADCZP2[2],corrADCZP2[3],corrADCZP2[4]);
588   printf(" ADCZEM1 [%1.2f] ADCZEM2 [%1.2f] \n",corrADCZEM1[0],corrADCZEM2[0]);
589   printf("*************************************************\n");*/
590     
591   // ---- Setting reco flags
592   UInt_t recoFlag=0;
593   UInt_t rFlags[32];
594   for(Int_t ifl=0; ifl<32; ifl++) rFlags[ifl]=0;
595
596   if(corrADCZN2[0]>fSignalThreshold)  rFlags[0] = 0x1;
597   if(corrADCZP2[0]>fSignalThreshold)  rFlags[1] = 0x1;
598   if(corrADCZEM1[0]>fSignalThreshold) rFlags[2] = 0x1;
599   if(corrADCZEM2[0]>fSignalThreshold) rFlags[3] = 0x1;
600   if(corrADCZN1[0]>fSignalThreshold)  rFlags[4] = 0x1;
601   if(corrADCZP1[0]>fSignalThreshold)  rFlags[5] = 0x1;
602   //
603   if(channelsOff==kTRUE)   rFlags[8] = 0x1;
604   if(chUnderflow == kTRUE) rFlags[9] = 0x1;
605   if(chOverflow==kTRUE)    rFlags[10] = 0x1;
606   recoFlag = rFlags[10] << 10 | rFlags[9] << 9 | rFlags[8] << 8 |
607              rFlags[5] << 5 | rFlags[4] << 4 | rFlags[3] << 3 |
608              rFlags[2] << 2 | rFlags[1] << 1 | rFlags[0];
609
610   // ******     Retrieving calibration data 
611   // --- Equalization coefficients ---------------------------------------------
612   Float_t equalCoeffZN1[5], equalCoeffZP1[5], equalCoeffZN2[5], equalCoeffZP2[5];
613   for(Int_t ji=0; ji<5; ji++){
614      equalCoeffZN1[ji] = fTowCalibData->GetZN1EqualCoeff(ji);
615      equalCoeffZP1[ji] = fTowCalibData->GetZP1EqualCoeff(ji); 
616      equalCoeffZN2[ji] = fTowCalibData->GetZN2EqualCoeff(ji); 
617      equalCoeffZP2[ji] = fTowCalibData->GetZP2EqualCoeff(ji); 
618   }
619   // --- Energy calibration factors ------------------------------------
620   Float_t calibEne[4];
621   // **** Energy calibration coefficient set to 1 
622   // **** (no trivial way to calibrate in p-p runs)
623   for(Int_t ij=0; ij<4; ij++) calibEne[ij] = fEnCalibData->GetEnCalib(ij);
624   
625   // ******     Equalization of detector responses
626   Float_t equalTowZN1[10], equalTowZN2[10], equalTowZP1[10], equalTowZP2[10];
627   for(Int_t gi=0; gi<10; gi++){
628      if(gi<5){
629        equalTowZN1[gi] = corrADCZN1[gi]*equalCoeffZN1[gi];
630        equalTowZP1[gi] = corrADCZP1[gi]*equalCoeffZP1[gi];
631        equalTowZN2[gi] = corrADCZN2[gi]*equalCoeffZN2[gi];
632        equalTowZP2[gi] = corrADCZP2[gi]*equalCoeffZP2[gi];
633      }
634      else{
635        equalTowZN1[gi] = corrADCZN1[gi]*equalCoeffZN1[gi-5];
636        equalTowZP1[gi] = corrADCZP1[gi]*equalCoeffZP1[gi-5];
637        equalTowZN2[gi] = corrADCZN2[gi]*equalCoeffZN2[gi-5];
638        equalTowZP2[gi] = corrADCZP2[gi]*equalCoeffZP2[gi-5];
639      }
640   }
641   // Ch. debug
642   /*printf("\n ------------- EQUALIZATION -------------\n");
643   printf(" ADCZN1 [%1.2f %1.2f %1.2f %1.2f %1.2f]\n",
644         equalTowZN1[0],equalTowZN1[1],equalTowZN1[2],equalTowZN1[3],equalTowZN1[4]);
645   printf(" ADCZP1 [%1.2f %1.2f %1.2f %1.2f %1.2f]\n",
646         equalTowZP1[0],equalTowZP1[1],equalTowZP1[2],equalTowZP1[3],equalTowZP1[4]);
647   printf(" ADCZN2 [%1.2f %1.2f %1.2f %1.2f %1.2f]\n",
648         equalTowZN2[0],equalTowZN2[1],equalTowZN2[2],equalTowZN2[3],equalTowZN2[4]);
649   printf(" ADCZP2 [%1.2f %1.2f %1.2f %1.2f %1.2f]\n",
650         equalTowZP2[0],equalTowZP2[1],equalTowZP2[2],equalTowZP2[3],equalTowZP2[4]);
651   printf(" ----------------------------------------\n");*/
652   
653   // ******     Summed response for hadronic calorimeter (SUMMED and then CALIBRATED!)
654   Float_t calibSumZN1[]={0,0}, calibSumZN2[]={0,0}, calibSumZP1[]={0,0}, calibSumZP2[]={0,0};
655   for(Int_t gi=0; gi<5; gi++){
656        calibSumZN1[0] += equalTowZN1[gi];
657        calibSumZP1[0] += equalTowZP1[gi];
658        calibSumZN2[0] += equalTowZN2[gi];
659        calibSumZP2[0] += equalTowZP2[gi];
660        //
661        calibSumZN1[1] += equalTowZN1[gi+5];
662        calibSumZP1[1] += equalTowZP1[gi+5];
663        calibSumZN2[1] += equalTowZN2[gi+5];
664        calibSumZP2[1] += equalTowZP2[gi+5];
665   }
666   // High gain chain
667   calibSumZN1[0] = calibSumZN1[0]*calibEne[0];
668   calibSumZP1[0] = calibSumZP1[0]*calibEne[1];
669   calibSumZN2[0] = calibSumZN2[0]*calibEne[2];
670   calibSumZP2[0] = calibSumZP2[0]*calibEne[3];
671   // Low gain chain
672   calibSumZN1[1] = calibSumZN1[1]*calibEne[0];
673   calibSumZP1[1] = calibSumZP1[1]*calibEne[1];
674   calibSumZN2[1] = calibSumZN2[1]*calibEne[2];
675   calibSumZP2[1] = calibSumZP2[1]*calibEne[3];
676   
677   // ******     Energy calibration of detector responses
678   Float_t calibTowZN1[10], calibTowZN2[10], calibTowZP1[10], calibTowZP2[10];
679   for(Int_t gi=0; gi<5; gi++){
680      // High gain chain
681      calibTowZN1[gi] = equalTowZN1[gi]*calibEne[0];
682      calibTowZP1[gi] = equalTowZP1[gi]*calibEne[1];
683      calibTowZN2[gi] = equalTowZN2[gi]*calibEne[2];
684      calibTowZP2[gi] = equalTowZP2[gi]*calibEne[3];
685      // Low gain chain
686      calibTowZN1[gi+5] = equalTowZN1[gi+5]*calibEne[0];
687      calibTowZP1[gi+5] = equalTowZP1[gi+5]*calibEne[1];
688      calibTowZN2[gi+5] = equalTowZN2[gi+5]*calibEne[2];
689      calibTowZP2[gi+5] = equalTowZP2[gi+5]*calibEne[3];
690   }
691   //
692   Float_t sumZEM[]={0,0}, calibZEM1[]={0,0}, calibZEM2[]={0,0};
693   calibZEM1[0] = corrADCZEM1[0]*calibEne[5];
694   calibZEM1[1] = corrADCZEM1[1]*calibEne[5];
695   calibZEM2[0] = corrADCZEM2[0]*calibEne[5];
696   calibZEM2[1] = corrADCZEM2[1]*calibEne[5];
697   for(Int_t k=0; k<2; k++) sumZEM[k] = calibZEM1[k] + calibZEM2[k];
698   // Ch. debug
699   /*printf("\n ------------- CALIBRATION -------------\n");
700   printf(" ADCZN1 [%1.2f %1.2f %1.2f %1.2f %1.2f]\n",
701         calibTowZN1[0],calibTowZN1[1],calibTowZN1[2],calibTowZN1[3],calibTowZN1[4]);
702   printf(" ADCZP1 [%1.2f %1.2f %1.2f %1.2f %1.2f]\n",
703         calibTowZP1[0],calibTowZP1[1],calibTowZP1[2],calibTowZP1[3],calibTowZP1[4]);
704   printf(" ADCZN2 [%1.2f %1.2f %1.2f %1.2f %1.2f]\n",
705         calibTowZN2[0],calibTowZN2[1],calibTowZN2[2],calibTowZN2[3],calibTowZN2[4]);
706   printf(" ADCZP2 [%1.2f %1.2f %1.2f %1.2f %1.2f]\n",
707         calibTowZP2[0],calibTowZP2[1],calibTowZP2[2],calibTowZP2[3],calibTowZP2[4]);
708   printf(" ADCZEM1 [%1.2f] ADCZEM2 [%1.2f] \n",calibZEM1[0],calibZEM2[0]);
709   printf(" ----------------------------------------\n");*/
710   
711   //  ******    No. of spectator and participants nucleons
712   //  Variables calculated to comply with ESD structure
713   //  *** N.B. -> They have a meaning only in Pb-Pb!!!!!!!!!!!!
714   Int_t nDetSpecNLeft=0, nDetSpecPLeft=0, nDetSpecNRight=0, nDetSpecPRight=0;
715   Int_t nGenSpec=0, nGenSpecLeft=0, nGenSpecRight=0;
716   Int_t nPart=0, nPartTotLeft=0, nPartTotRight=0;
717   Double_t impPar=0., impPar1=0., impPar2=0.;
718   
719   // create the output tree
720   AliZDCReco* reco = new AliZDCReco(calibSumZN1, calibSumZP1, calibSumZN2, calibSumZP2, 
721                    calibTowZN1, calibTowZP1, calibTowZN2, calibTowZP2, 
722                    calibZEM1, calibZEM2, sPMRef1, sPMRef2,
723                    nDetSpecNLeft, nDetSpecPLeft, nDetSpecNRight, nDetSpecPRight, 
724                    nGenSpec, nGenSpecLeft, nGenSpecRight, 
725                    nPart, nPartTotLeft, nPartTotRight, 
726                    impPar, impPar1, impPar2,
727                    recoFlag);
728                   
729   const Int_t kBufferSize = 4000;
730   clustersTree->Branch("ZDC", "AliZDCReco", &reco, kBufferSize);
731   // write the output tree
732   clustersTree->Fill();
733 }
734
735 //_____________________________________________________________________________
736 void AliZDCReconstructor::ReconstructEventPbPb(TTree *clustersTree, 
737         Float_t* corrADCZN1, Float_t* corrADCZP1, Float_t* corrADCZN2, Float_t* corrADCZP2,
738         Float_t* corrADCZEM1, Float_t* corrADCZEM2, Float_t* sPMRef1, Float_t* sPMRef2,
739         Bool_t channelsOff, Bool_t chUnderflow, Bool_t chOverflow) const 
740 {
741   // ****************** Reconstruct one event ******************
742   
743   // ---- Setting reco flags
744   UInt_t recoFlag=0;
745   UInt_t rFlags[32];
746   for(Int_t ifl=0; ifl<32; ifl++) rFlags[ifl]=0;
747
748   if(corrADCZN2[0]>fSignalThreshold)  rFlags[0] = 0x1;
749   if(corrADCZP2[0]>fSignalThreshold)  rFlags[1] = 0x1;
750   if(corrADCZEM1[0]>fSignalThreshold) rFlags[2] = 0x1;
751   if(corrADCZEM2[0]>fSignalThreshold) rFlags[3] = 0x1;
752   if(corrADCZN1[0]>fSignalThreshold)  rFlags[4] = 0x1;
753   if(corrADCZP1[0]>fSignalThreshold)  rFlags[5] = 0x1;
754   //
755   if(channelsOff==kTRUE)   rFlags[8] = 0x1;
756   if(chUnderflow == kTRUE) rFlags[9] = 0x1;
757   if(chOverflow==kTRUE)    rFlags[10] = 0x1;
758   recoFlag = rFlags[10] << 10 | rFlags[9] << 9 | rFlags[8] << 8 |
759              rFlags[5] << 5 | rFlags[4] << 4 | rFlags[3] << 3 |
760              rFlags[2] << 2 | rFlags[1] << 1 | rFlags[0];
761
762   // ******     Retrieving calibration data 
763   // --- Equalization coefficients ---------------------------------------------
764   Float_t equalCoeffZN1[5], equalCoeffZP1[5], equalCoeffZN2[5], equalCoeffZP2[5];
765   for(Int_t ji=0; ji<5; ji++){
766      equalCoeffZN1[ji] = fTowCalibData->GetZN1EqualCoeff(ji);
767      equalCoeffZP1[ji] = fTowCalibData->GetZP1EqualCoeff(ji); 
768      equalCoeffZN2[ji] = fTowCalibData->GetZN2EqualCoeff(ji); 
769      equalCoeffZP2[ji] = fTowCalibData->GetZP2EqualCoeff(ji); 
770   }
771   // --- Energy calibration factors ------------------------------------
772   Float_t valFromOCDB[6], calibEne[6];
773   for(Int_t ij=0; ij<6; ij++){
774     valFromOCDB[ij] = fEnCalibData->GetEnCalib(ij);
775     if(ij<4){
776       if(valFromOCDB[ij]!=0) calibEne[ij] = fBeamEnergy/valFromOCDB[ij];
777       else AliWarning(" Value from OCDB for E calibration = 0 !!!\n");
778     } 
779     else calibEne[ij] = valFromOCDB[ij];
780   }
781   
782   // ******     Equalization of detector responses
783   Float_t equalTowZN1[10], equalTowZN2[10], equalTowZP1[10], equalTowZP2[10];
784   for(Int_t gi=0; gi<10; gi++){
785      if(gi<5){
786        equalTowZN1[gi] = corrADCZN1[gi]*equalCoeffZN1[gi];
787        equalTowZP1[gi] = corrADCZP1[gi]*equalCoeffZP1[gi];
788        equalTowZN2[gi] = corrADCZN2[gi]*equalCoeffZN2[gi];
789        equalTowZP2[gi] = corrADCZP2[gi]*equalCoeffZP2[gi];
790      }
791      else{
792        equalTowZN1[gi] = corrADCZN1[gi]*equalCoeffZN1[gi-5];
793        equalTowZP1[gi] = corrADCZP1[gi]*equalCoeffZP1[gi-5];
794        equalTowZN2[gi] = corrADCZN2[gi]*equalCoeffZN2[gi-5];
795        equalTowZP2[gi] = corrADCZP2[gi]*equalCoeffZP2[gi-5];
796      }
797   }
798   
799   // ******     Summed response for hadronic calorimeter (SUMMED and then CALIBRATED!)
800   Float_t calibSumZN1[]={0,0}, calibSumZN2[]={0,0}, calibSumZP1[]={0,0}, calibSumZP2[]={0,0};
801   for(Int_t gi=0; gi<5; gi++){
802        calibSumZN1[0] += equalTowZN1[gi];
803        calibSumZP1[0] += equalTowZP1[gi];
804        calibSumZN2[0] += equalTowZN2[gi];
805        calibSumZP2[0] += equalTowZP2[gi];
806        //
807        calibSumZN1[1] += equalTowZN1[gi+5];
808        calibSumZP1[1] += equalTowZP1[gi+5];
809        calibSumZN2[1] += equalTowZN2[gi+5];
810        calibSumZP2[1] += equalTowZP2[gi+5];
811   }
812   // High gain chain
813   calibSumZN1[0] = calibSumZN1[0]*calibEne[0]/8.;
814   calibSumZP1[0] = calibSumZP1[0]*calibEne[1]/8.;
815   calibSumZN2[0] = calibSumZN2[0]*calibEne[2]/8.;
816   calibSumZP2[0] = calibSumZP2[0]*calibEne[3]/8.;
817   // Low gain chain
818   calibSumZN1[1] = calibSumZN1[1]*calibEne[0];
819   calibSumZP1[1] = calibSumZP1[1]*calibEne[1];
820   calibSumZN2[1] = calibSumZN2[1]*calibEne[2];
821   calibSumZP2[1] = calibSumZP2[1]*calibEne[3];
822   //
823   Float_t sumZEM[]={0,0}, calibZEM1[]={0,0}, calibZEM2[]={0,0};
824   calibZEM1[0] = corrADCZEM1[0]*calibEne[5]/8.;
825   calibZEM1[1] = corrADCZEM1[1]*calibEne[5];
826   calibZEM2[0] = corrADCZEM2[0]*calibEne[5]/8.;
827   calibZEM2[1] = corrADCZEM2[1]*calibEne[5];
828   for(Int_t k=0; k<2; k++) sumZEM[k] = calibZEM1[k] + calibZEM2[k];
829     
830   // ******     Energy calibration of detector responses
831   Float_t calibTowZN1[10], calibTowZN2[10], calibTowZP1[10], calibTowZP2[10];
832   for(Int_t gi=0; gi<5; gi++){
833      // High gain chain
834      calibTowZN1[gi] = equalTowZN1[gi]*calibEne[0]/8.;
835      calibTowZP1[gi] = equalTowZP1[gi]*calibEne[1]/8.;
836      calibTowZN2[gi] = equalTowZN2[gi]*calibEne[2]/8.;
837      calibTowZP2[gi] = equalTowZP2[gi]*calibEne[3]/8.;
838      // Low gain chain
839      calibTowZN1[gi+5] = equalTowZN1[gi+5]*calibEne[0];
840      calibTowZP1[gi+5] = equalTowZP1[gi+5]*calibEne[1];
841      calibTowZN2[gi+5] = equalTowZN2[gi+5]*calibEne[2];
842      calibTowZP2[gi+5] = equalTowZP2[gi+5]*calibEne[3];
843   }
844   
845   //  ******    Number of detected spectator nucleons
846   Int_t nDetSpecNLeft=0, nDetSpecPLeft=0, nDetSpecNRight=0, nDetSpecPRight=0;
847   if(fBeamEnergy!=0){
848     nDetSpecNLeft = (Int_t) (calibSumZN1[0]/fBeamEnergy);
849     nDetSpecPLeft = (Int_t) (calibSumZP1[0]/fBeamEnergy);
850     nDetSpecNRight = (Int_t) (calibSumZN2[0]/fBeamEnergy);
851     nDetSpecPRight = (Int_t) (calibSumZP2[0]/fBeamEnergy);
852   }
853   else AliWarning(" ATTENTION!!! fBeamEnergy=0 -> N_spec will be ZERO!!! \n");
854   /*printf("\n\t AliZDCReconstructor -> nDetSpecNLeft %d, nDetSpecPLeft %d,"
855     " nDetSpecNRight %d, nDetSpecPRight %d\n",nDetSpecNLeft, nDetSpecPLeft, 
856     nDetSpecNRight, nDetSpecPRight);*/
857   
858   Int_t nGenSpec=0, nGenSpecA=0, nGenSpecC=0;
859   Int_t nPart=0, nPartA=0, nPartC=0;
860   Double_t b=0., bA=0., bC=0.;
861   
862   if(fIsCalibrationMB == kFALSE){
863     // ******   Reconstruction parameters ------------------ 
864     // Ch. debug
865     //fRecoParam->Print("");
866     //
867  
868     if (!fRecoParam) fRecoParam = const_cast<AliZDCRecoParam*>(GetRecoParam()); 
869  
870     TH2F *hZDCvsZEM  = fRecoParam->GethZDCvsZEM();
871     TH2F *hZDCCvsZEM = fRecoParam->GethZDCCvsZEM();
872     TH2F *hZDCAvsZEM = fRecoParam->GethZDCAvsZEM();
873     TH1D *hNpartDist = fRecoParam->GethNpartDist();
874     TH1D *hbDist = fRecoParam->GethbDist();    
875     Float_t ClkCenter = fRecoParam->GetClkCenter();
876     //
877     Double_t xHighEdge = hZDCvsZEM->GetXaxis()->GetXmax();
878     Double_t origin = xHighEdge*ClkCenter;
879     // Ch. debug
880     //printf("\n\n  xHighEdge %1.2f, origin %1.4f \n", xHighEdge, origin);
881     //
882     // ====> Summed ZDC info (sideA+side C)
883     TF1 *line = new TF1("line","[0]*x+[1]",0.,xHighEdge);
884     Float_t y = (calibSumZN1[0]+calibSumZP1[0]+calibSumZN2[0]+calibSumZP2[0])/1000.;
885     Float_t x = (calibZEM1[0]+calibZEM2[0])/1000.;
886     line->SetParameter(0, y/(x-origin));
887     line->SetParameter(1, -origin*y/(x-origin));
888     // Ch. debug
889     //printf("  ***************** Summed ZDC info (sideA+side C) \n");
890     //printf("  E_{ZEM} %1.4f, E_{ZDC} %1.2f, TF1: %1.2f*x + %1.2f   ", x, y,y/(x-origin),-origin*y/(x-origin));
891     //
892     Double_t countPerc=0;
893     Double_t xBinCenter=0, yBinCenter=0;
894     for(Int_t nbinx=1; nbinx<=hZDCvsZEM->GetNbinsX(); nbinx++){
895       for(Int_t nbiny=1; nbiny<=hZDCvsZEM->GetNbinsY(); nbiny++){
896          xBinCenter = hZDCvsZEM->GetXaxis()->GetBinCenter(nbinx);
897          yBinCenter = hZDCvsZEM->GetYaxis()->GetBinCenter(nbiny);
898          //
899          if(line->GetParameter(0)>0){
900            if(yBinCenter < (line->GetParameter(0)*xBinCenter + line->GetParameter(1))){
901              countPerc += hZDCvsZEM->GetBinContent(nbinx,nbiny);
902              // Ch. debug
903              /*printf(" xBinCenter  %1.3f, yBinCenter %1.0f,  countPerc %1.0f\n", 
904                 xBinCenter, yBinCenter, countPerc);*/
905            }
906          }
907          else{
908            if(yBinCenter > (line->GetParameter(0)*xBinCenter + line->GetParameter(1))){
909              countPerc += hZDCvsZEM->GetBinContent(nbinx,nbiny);
910              // Ch. debug
911              /*printf(" xBinCenter  %1.3f, yBinCenter %1.0f,  countPerc %1.0f\n", 
912                 xBinCenter, yBinCenter, countPerc);*/
913            }
914          }
915       }
916     }
917     //
918     Double_t xSecPerc = 0.;
919     if(hZDCvsZEM->GetEntries()!=0){ 
920       xSecPerc = countPerc/hZDCvsZEM->GetEntries();
921     }
922     else{
923       AliWarning("  Histogram hZDCvsZEM from OCDB has no entries!!!");
924     }
925     // Ch. debug
926     //printf("  xSecPerc %1.4f  \n", xSecPerc);
927
928     // ====> side C
929     TF1 *lineC = new TF1("lineC","[0]*x+[1]",0.,xHighEdge);
930     Float_t yC = (calibSumZN1[0]+calibSumZP1[0])/1000.;
931     lineC->SetParameter(0, yC/(x-origin));
932     lineC->SetParameter(1, -origin*yC/(x-origin));
933     // Ch. debug
934     //printf("  ***************** Side C \n");
935     //printf("  E_{ZEM} %1.4f, E_{ZDCC} %1.2f, TF1: %1.2f*x + %1.2f   ", x, yC,yC/(x-origin),-origin*yC/(x-origin));
936     //
937     Double_t countPercC=0;
938     Double_t xBinCenterC=0, yBinCenterC=0;
939     for(Int_t nbinx=1; nbinx<=hZDCCvsZEM->GetNbinsX(); nbinx++){
940       for(Int_t nbiny=1; nbiny<=hZDCCvsZEM->GetNbinsY(); nbiny++){
941          xBinCenterC = hZDCCvsZEM->GetXaxis()->GetBinCenter(nbinx);
942          yBinCenterC = hZDCCvsZEM->GetYaxis()->GetBinCenter(nbiny);
943          if(lineC->GetParameter(0)>0){
944            if(yBinCenterC < (lineC->GetParameter(0)*xBinCenterC + lineC->GetParameter(1))){
945              countPercC += hZDCCvsZEM->GetBinContent(nbinx,nbiny);
946            }
947          }
948          else{
949            if(yBinCenterC > (lineC->GetParameter(0)*xBinCenterC + lineC->GetParameter(1))){
950              countPercC += hZDCCvsZEM->GetBinContent(nbinx,nbiny);
951            }
952          }
953       }
954     }
955     //
956     Double_t xSecPercC = 0.;
957     if(hZDCCvsZEM->GetEntries()!=0){ 
958       xSecPercC = countPercC/hZDCCvsZEM->GetEntries();
959     }
960     else{
961       AliWarning("  Histogram hZDCCvsZEM from OCDB has no entries!!!");
962     }
963     // Ch. debug
964     //printf("  xSecPercC %1.4f  \n", xSecPercC);
965     
966     // ====> side A
967     TF1 *lineA = new TF1("lineA","[0]*x+[1]",0.,xHighEdge);
968     Float_t yA = (calibSumZN2[0]+calibSumZP2[0])/1000.;
969     lineA->SetParameter(0, yA/(x-origin));
970     lineA->SetParameter(1, -origin*yA/(x-origin));
971     //
972     // Ch. debug
973     //printf("  ***************** Side A \n");
974     //printf("  E_{ZEM} %1.4f, E_{ZDCA} %1.2f, TF1: %1.2f*x + %1.2f   ", x, yA,yA/(x-origin),-origin*yA/(x-origin));
975     //
976     Double_t countPercA=0;
977     Double_t xBinCenterA=0, yBinCenterA=0;
978     for(Int_t nbinx=1; nbinx<=hZDCAvsZEM->GetNbinsX(); nbinx++){
979       for(Int_t nbiny=1; nbiny<=hZDCAvsZEM->GetNbinsY(); nbiny++){
980          xBinCenterA = hZDCAvsZEM->GetXaxis()->GetBinCenter(nbinx);
981          yBinCenterA = hZDCAvsZEM->GetYaxis()->GetBinCenter(nbiny);
982          if(lineA->GetParameter(0)>0){
983            if(yBinCenterA < (lineA->GetParameter(0)*xBinCenterA + lineA->GetParameter(1))){
984              countPercA += hZDCAvsZEM->GetBinContent(nbinx,nbiny);
985            }
986          }
987          else{
988            if(yBinCenterA > (lineA->GetParameter(0)*xBinCenterA + lineA->GetParameter(1))){
989              countPercA += hZDCAvsZEM->GetBinContent(nbinx,nbiny);
990            }
991          }
992       }
993     }
994     //
995     Double_t xSecPercA = 0.;
996     if(hZDCAvsZEM->GetEntries()!=0){ 
997       xSecPercA = countPercA/hZDCAvsZEM->GetEntries();
998     }
999     else{
1000       AliWarning("  Histogram hZDCAvsZEM from OCDB has no entries!!!");
1001     }
1002     // Ch. debug
1003     //printf("  xSecPercA %1.4f  \n", xSecPercA);
1004     
1005     //  ******    Number of participants (from E_ZDC vs. E_ZEM correlation)
1006     Int_t nPart=0, nPartC=0, nPartA=0;
1007     Double_t nPartFrac=0., nPartFracC=0., nPartFracA=0.;
1008     for(Int_t npbin=1; npbin<hNpartDist->GetNbinsX(); npbin++){
1009       nPartFrac += (hNpartDist->GetBinContent(npbin))/(hNpartDist->GetEntries());
1010       if((1.-nPartFrac) < xSecPerc){
1011         nPart = (Int_t) hNpartDist->GetBinLowEdge(npbin);
1012         // Ch. debug
1013         //printf("  ***************** Summed ZDC info (sideA+side C) \n");
1014         //printf("  nPartFrac %1.4f, nPart %d\n", nPartFrac, nPart);
1015         break;
1016       }
1017     }
1018     if(nPart<0) nPart=0;
1019     //
1020     for(Int_t npbin=1; npbin<hNpartDist->GetNbinsX(); npbin++){
1021       nPartFracC += (hNpartDist->GetBinContent(npbin))/(hNpartDist->GetEntries());
1022       if((1.-nPartFracC) < xSecPercC){
1023         nPartC = (Int_t) hNpartDist->GetBinLowEdge(npbin);
1024         // Ch. debug
1025         //printf("  ***************** Side C \n");
1026         //printf("  nPartFracC %1.4f, nPartC %d\n", nPartFracC, nPartC);
1027         break;
1028     }
1029     }
1030     if(nPartC<0) nPartC=0;
1031     //
1032     for(Int_t npbin=1; npbin<hNpartDist->GetNbinsX(); npbin++){
1033       nPartFracA += (hNpartDist->GetBinContent(npbin))/(hNpartDist->GetEntries());
1034       if((1.-nPartFracA) < xSecPercA){
1035         nPartA = (Int_t) hNpartDist->GetBinLowEdge(npbin);
1036         // Ch. debug
1037         //printf("  ***************** Side A \n");
1038         //printf("  nPartFracA %1.4f, nPartA %d\n\n", nPartFracA, nPartA);
1039         break;
1040       }
1041     }
1042     if(nPartA<0) nPartA=0;
1043     
1044     //  ******    Impact parameter (from E_ZDC vs. E_ZEM correlation)
1045     Float_t b=0, bC=0, bA=0;
1046     Double_t bFrac=0., bFracC=0., bFracA=0.;
1047     for(Int_t ibbin=1; ibbin<hbDist->GetNbinsX(); ibbin++){
1048       bFrac += (hbDist->GetBinContent(ibbin))/(hbDist->GetEntries());
1049       if((1.-bFrac) < xSecPerc){
1050         b = hbDist->GetBinLowEdge(ibbin);
1051         break;
1052       }
1053     }
1054     //
1055     for(Int_t ibbin=1; ibbin<hbDist->GetNbinsX(); ibbin++){
1056       bFracC += (hbDist->GetBinContent(ibbin))/(hbDist->GetEntries());
1057       if((1.-bFracC) < xSecPercC){
1058         bC = hbDist->GetBinLowEdge(ibbin);
1059         break;
1060       }
1061     }
1062     //
1063     for(Int_t ibbin=1; ibbin<hbDist->GetNbinsX(); ibbin++){
1064       bFracA += (hbDist->GetBinContent(ibbin))/(hNpartDist->GetEntries());
1065       if((1.-bFracA) < xSecPercA){
1066         bA = hbDist->GetBinLowEdge(ibbin);
1067         break;
1068       }
1069     }
1070
1071     //  ******  Number of spectator nucleons 
1072     Int_t nGenSpec=0, nGenSpecC=0, nGenSpecA=0;
1073     //
1074     nGenSpec = 416 - nPart;
1075     nGenSpecC = 416 - nPartC;
1076     nGenSpecA = 416 - nPartA;
1077     if(nGenSpec>416) nGenSpec=416; if(nGenSpec<0) nGenSpec=0;
1078     if(nGenSpecC>416) nGenSpecC=416; if(nGenSpecC<0) nGenSpecC=0;
1079     if(nGenSpecA>416) nGenSpecA=416; if(nGenSpecA<0) nGenSpecA=0;
1080   
1081     //  Ch. debug
1082     /*printf("\n\t AliZDCReconstructor -> calibSumZN1[0] %1.0f, calibSumZP1[0] %1.0f,"
1083         "  calibSumZN2[0] %1.0f, calibSumZP2[0] %1.0f, corrADCZEMHG %1.0f\n", 
1084         calibSumZN1[0],calibSumZP1[0],calibSumZN2[0],calibSumZP2[0],corrADCZEMHG);
1085     printf("\t AliZDCReconstructor -> nGenSpecLeft %d nGenSpecRight %d\n", 
1086         nGenSpecLeft, nGenSpecRight);
1087     printf("\t AliZDCReconstructor ->  NpartL %d,  NpartR %d,  b %1.2f fm\n\n",nPartTotLeft, nPartTotRight, impPar);
1088     */
1089   
1090     delete lineC;  delete lineA;
1091
1092   } // ONLY IF fIsCalibrationMB==kFALSE
1093   
1094   AliZDCReco* reco = new AliZDCReco(calibSumZN1, calibSumZP1, calibSumZN2, calibSumZP2, 
1095                   calibTowZN1, calibTowZP1, calibTowZN2, calibTowZP2, 
1096                   calibZEM1, calibZEM2, sPMRef1, sPMRef2,
1097                   nDetSpecNLeft, nDetSpecPLeft, nDetSpecNRight, nDetSpecPRight, 
1098                   nGenSpec, nGenSpecA, nGenSpecC, 
1099                   nPart, nPartA, nPartC, b, bA, bC,
1100                   recoFlag);
1101                     
1102   const Int_t kBufferSize = 4000;
1103   clustersTree->Branch("ZDC", "AliZDCReco", &reco, kBufferSize);
1104   reco->Print("");
1105   // write the output tree
1106   clustersTree->Fill();
1107 }
1108
1109 //_____________________________________________________________________________
1110 void AliZDCReconstructor::BuildRecoParam(Float_t ZDCC, Float_t ZDCA, Float_t ZEM) const
1111 {
1112   // Calculate RecoParam object for Pb-Pb data
1113   (fRecoParam->GethZDCvsZEM())->Fill(ZDCC+ZDCA, ZEM);
1114   (fRecoParam->GethZDCCvsZEM())->Fill(ZDCC, ZEM);
1115   (fRecoParam->GethZDCAvsZEM())->Fill(ZDCA, ZEM);
1116  
1117 }
1118
1119 //_____________________________________________________________________________
1120 void AliZDCReconstructor::FillZDCintoESD(TTree *clustersTree, AliESDEvent* esd) const
1121 {
1122   // fill energies and number of participants to the ESD
1123
1124   if(fIsCalibrationMB==kTRUE) WritePbPbRecoParamInOCDB();
1125   
1126   AliZDCReco reco;
1127   AliZDCReco* preco = &reco;
1128   clustersTree->SetBranchAddress("ZDC", &preco);
1129
1130   clustersTree->GetEntry(0);
1131   //
1132   AliESDZDC * esdzdc = esd->GetESDZDC();
1133   Float_t tZN1Ene[5], tZN2Ene[5], tZP1Ene[5], tZP2Ene[5];
1134   Float_t tZN1EneLR[5], tZN2EneLR[5], tZP1EneLR[5], tZP2EneLR[5];
1135   for(Int_t i=0; i<5; i++){
1136      tZN1Ene[i] = reco.GetZN1HREnTow(i);
1137      tZN2Ene[i] = reco.GetZN2HREnTow(i);
1138      tZP1Ene[i] = reco.GetZP1HREnTow(i);
1139      tZP2Ene[i] = reco.GetZP2HREnTow(i);
1140      //
1141      tZN1EneLR[i] = reco.GetZN1LREnTow(i);
1142      tZN2EneLR[i] = reco.GetZN2LREnTow(i);
1143      tZP1EneLR[i] = reco.GetZP1LREnTow(i);
1144      tZP2EneLR[i] = reco.GetZP2LREnTow(i);
1145   }
1146   //
1147   esdzdc->SetZN1TowerEnergy(tZN1Ene);
1148   esdzdc->SetZN2TowerEnergy(tZN2Ene);
1149   esdzdc->SetZP1TowerEnergy(tZP1Ene);
1150   esdzdc->SetZP2TowerEnergy(tZP2Ene);
1151   //
1152   esdzdc->SetZN1TowerEnergyLR(tZN1EneLR);
1153   esdzdc->SetZN2TowerEnergyLR(tZN2EneLR);
1154   esdzdc->SetZP1TowerEnergyLR(tZP1EneLR);
1155   esdzdc->SetZP2TowerEnergyLR(tZP2EneLR);
1156   // 
1157   Int_t nPart  = reco.GetNParticipants();
1158   Int_t nPartA = reco.GetNPartSideA();
1159   Int_t nPartC = reco.GetNPartSideC();
1160   Double_t b  = reco.GetImpParameter();
1161   Double_t bA = reco.GetImpParSideA();
1162   Double_t bC = reco.GetImpParSideC();
1163   UInt_t recoFlag = reco.GetRecoFlag();
1164   //
1165   esd->SetZDC(reco.GetZN1HREnergy(), reco.GetZP1HREnergy(), reco.GetZEM1HRsignal(), 
1166               reco.GetZEM2HRsignal(), reco.GetZN2HREnergy(), reco.GetZP2HREnergy(), 
1167               nPart, nPartA, nPartC, b, bA, bC, recoFlag);
1168   
1169 }
1170
1171 //_____________________________________________________________________________
1172 AliCDBStorage* AliZDCReconstructor::SetStorage(const char *uri) 
1173 {
1174   // Setting the storage
1175
1176   Bool_t deleteManager = kFALSE;
1177   
1178   AliCDBManager *manager = AliCDBManager::Instance();
1179   AliCDBStorage *defstorage = manager->GetDefaultStorage();
1180   
1181   if(!defstorage || !(defstorage->Contains("ZDC"))){ 
1182      AliWarning("No default storage set or default storage doesn't contain ZDC!");
1183      manager->SetDefaultStorage(uri);
1184      deleteManager = kTRUE;
1185   }
1186  
1187   AliCDBStorage *storage = manager->GetDefaultStorage();
1188
1189   if(deleteManager){
1190     AliCDBManager::Instance()->UnsetDefaultStorage();
1191     defstorage = 0;   // the storage is killed by AliCDBManager::Instance()->Destroy()
1192   }
1193
1194   return storage; 
1195 }
1196
1197 //_____________________________________________________________________________
1198 AliZDCPedestals* AliZDCReconstructor::GetPedestalData() const
1199 {
1200
1201   // Getting pedestal calibration object for ZDC set
1202
1203   AliCDBEntry  *entry = AliCDBManager::Instance()->Get("ZDC/Calib/Pedestals");
1204   if(!entry) AliFatal("No calibration data loaded!");  
1205
1206   AliZDCPedestals *calibdata = dynamic_cast<AliZDCPedestals*>  (entry->GetObject());
1207   if(!calibdata)  AliFatal("Wrong calibration object in calibration  file!");
1208
1209   return calibdata;
1210 }
1211
1212 //_____________________________________________________________________________
1213 AliZDCEnCalib* AliZDCReconstructor::GetEnergyCalibData() const
1214 {
1215
1216   // Getting energy and equalization calibration object for ZDC set
1217
1218   AliCDBEntry  *entry = AliCDBManager::Instance()->Get("ZDC/Calib/EnergyCalib");
1219   if(!entry) AliFatal("No calibration data loaded!");  
1220
1221   AliZDCEnCalib *calibdata = dynamic_cast<AliZDCEnCalib*> (entry->GetObject());
1222   if(!calibdata)  AliFatal("Wrong calibration object in calibration  file!");
1223
1224   return calibdata;
1225 }
1226
1227 //_____________________________________________________________________________
1228 AliZDCTowerCalib* AliZDCReconstructor::GetTowerCalibData() const
1229 {
1230
1231   // Getting energy and equalization calibration object for ZDC set
1232
1233   AliCDBEntry  *entry = AliCDBManager::Instance()->Get("ZDC/Calib/TowerCalib");
1234   if(!entry) AliFatal("No calibration data loaded!");  
1235
1236   AliZDCTowerCalib *calibdata = dynamic_cast<AliZDCTowerCalib*> (entry->GetObject());
1237   if(!calibdata)  AliFatal("Wrong calibration object in calibration  file!");
1238
1239   return calibdata;
1240 }
1241
1242 //_____________________________________________________________________________
1243 AliZDCRecoParampp* AliZDCReconstructor::GetppRecoParamFromOCDB() const
1244 {
1245
1246   // Getting reconstruction parameters from OCDB
1247
1248   AliCDBEntry  *entry = AliCDBManager::Instance()->Get("ZDC/Calib/RecoParampp");
1249   if(!entry) AliFatal("No RecoParam data found in OCDB!");  
1250   
1251   AliZDCRecoParampp *param = dynamic_cast<AliZDCRecoParampp*> (entry->GetObject());
1252   if(!param)  AliFatal("No RecoParam object in OCDB entry!");
1253   
1254   return param;
1255
1256 }
1257
1258 //_____________________________________________________________________________
1259 AliZDCRecoParamPbPb* AliZDCReconstructor::GetPbPbRecoParamFromOCDB() const
1260 {
1261
1262   // Getting reconstruction parameters from OCDB
1263
1264   AliCDBEntry  *entry = AliCDBManager::Instance()->Get("ZDC/Calib/RecoParamPbPb");
1265   if(!entry) AliFatal("No RecoParam data found in OCDB!");  
1266   
1267   AliZDCRecoParamPbPb *param = dynamic_cast<AliZDCRecoParamPbPb*> (entry->GetObject());
1268   if(!param)  AliFatal("No RecoParam object in OCDB entry!");
1269   
1270   return param;
1271
1272 }
1273
1274 //_____________________________________________________________________________
1275 void AliZDCReconstructor::WritePbPbRecoParamInOCDB() const
1276 {
1277
1278   // Writing Pb-Pb reconstruction parameters from OCDB
1279
1280   AliCDBManager *man = AliCDBManager::Instance();
1281   AliCDBMetaData *md= new AliCDBMetaData();
1282   md->SetResponsible("Chiara Oppedisano");
1283   md->SetComment("ZDC Pb-Pb reconstruction parameters");
1284   md->SetObjectClassName("AliZDCRecoParamPbPb");
1285   AliCDBId id("ZDC/Calib/RecoParamPbPb",fNRun,AliCDBRunRange::Infinity());
1286   man->Put(fRecoParam, id, md);
1287
1288 }
1289