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