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