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