ZDC/AliZDCReconstructor.cxx

   1 /**************************************************************************
   2  * Copyright(c) 1998-1999, ALICE Experiment at CERN, All rights reserved. *
   3  *                                                                        *
   4  * Author: The ALICE Off-line Project.                                    *
   5  * Contributors are mentioned in the code where appropriate.              *
   6  *                                                                        *
   7  * Permission to use, copy, modify and distribute this software and its   *
   8  * documentation strictly for non-commercial purposes is hereby granted   *
   9  * without fee, provided that the above copyright notice appears in all   *
  10  * copies and that both the copyright notice and this permission notice   *
  11  * appear in the supporting documentation. The authors make no claims     *
  12  * about the suitability of this software for any purpose. It is          *
  13  * provided "as is" without express or implied warranty.                  *
  14  **************************************************************************/
  15
  16 /* $Id$ */
  17
  18 ///////////////////////////////////////////////////////////////////////////////
  19 //                                                                           //
  20 // class for ZDC reconstruction                                              //
  21 //                                                                           //
  22 ///////////////////////////////////////////////////////////////////////////////
  23
  24
  25 #include <TF1.h>
  26
  27 #include "AliRunLoader.h"
  28 #include "AliRawReader.h"
  29 #include "AliESDEvent.h"
  30 #include "AliZDCDigit.h"
  31 #include "AliZDCRawStream.h"
  32 #include "AliZDCReco.h"
  33 #include "AliZDCReconstructor.h"
  34 #include "AliZDCCalibData.h"
  35
  36
  37 ClassImp(AliZDCReconstructor)
  38
  39
  40 //_____________________________________________________________________________
  41 AliZDCReconstructor:: AliZDCReconstructor() :
  42
  43   fZNCen(new TF1("fZNCen",
  44         "(-2.287920+sqrt(2.287920*2.287920-4*(-0.007629)*(11.921710-x)))/(2*(-0.007629))",0.,164.)),
  45   fZNPer(new TF1("fZNPer",
  46       "(-37.812280-sqrt(37.812280*37.812280-4*(-0.190932)*(-1709.249672-x)))/(2*(-0.190932))",0.,164.)),
  47   fZPCen(new TF1("fZPCen",
  48        "(-1.321353+sqrt(1.321353*1.321353-4*(-0.007283)*(3.550697-x)))/(2*(-0.007283))",0.,60.)),
  49   fZPPer(new TF1("fZPPer",
  50       "(-42.643308-sqrt(42.643308*42.643308-4*(-0.310786)*(-1402.945615-x)))/(2*(-0.310786))",0.,60.)),
  51   fZDCCen(new TF1("fZDCCen",
  52       "(-1.934991+sqrt(1.934991*1.934991-4*(-0.004080)*(15.111124-x)))/(2*(-0.004080))",0.,225.)),
  53   fZDCPer(new TF1("fZDCPer",
  54       "(-34.380639-sqrt(34.380639*34.380639-4*(-0.104251)*(-2612.189017-x)))/(2*(-0.104251))",0.,225.)),
  55   fbCen(new TF1("fbCen","-0.056923+0.079703*x-0.0004301*x*x+0.000001366*x*x*x",0.,220.)),
  56   fbPer(new TF1("fbPer","17.943998-0.046846*x+0.000074*x*x",0.,220.)),
  57   fZEMn(new TF1("fZEMn","126.2-0.05399*x+0.000005679*x*x",0.,4000.)),
  58   fZEMp(new TF1("fZEMp","82.49-0.03611*x+0.00000385*x*x",0.,4000.)),
  59   fZEMsp(new TF1("fZEMsp","208.7-0.09006*x+0.000009526*x*x",0.,4000.)),
  60   fZEMb(new TF1("fZEMb",
  61         "16.06-0.01633*x+1.44e-5*x*x-6.778e-9*x*x*x+1.438e-12*x*x*x*x-1.112e-16*x*x*x*x*x",0.,4000.)),
  62   //
  63   fCalibData(GetCalibData())
  64
  65 {
  66   // **** Default constructor
  67
  68 }
  69
  70
  71 //_____________________________________________________________________________
  72 AliZDCReconstructor::~AliZDCReconstructor()
  73 {
  74 // destructor
  75
  76   delete fZNCen;
  77   delete fZNPer;
  78   delete fZPCen;
  79   delete fZPPer;
  80   delete fZDCCen;
  81   delete fZDCPer;
  82   delete fbCen;
  83   delete fbPer;
  84   delete fZEMn;
  85   delete fZEMp;
  86   delete fZEMsp;
  87   delete fZEMb;
  88
  89 }
  90
  91
  92 //_____________________________________________________________________________
  93 void AliZDCReconstructor::Reconstruct(TTree* digitsTree, TTree* clustersTree) const
  94 {
  95   // *** Local ZDC reconstruction for digits
  96   // Works on the current event
  97
  98   // Retrieving calibration data
  99   Float_t meanPed[47];
 100   for(Int_t jj=0; jj<47; jj++) meanPed[jj] = fCalibData->GetMeanPed(jj);
 101
 102   // get digits
 103   AliZDCDigit digit;
 104   AliZDCDigit* pdigit = &digit;
 105   digitsTree->SetBranchAddress("ZDC", &pdigit);
 106
 107   // loop over digits
 108     Float_t tZN1CorrHG[]={0.,0.,0.,0.,0.}, tZP1CorrHG[]={0.,0.,0.,0.,0.};
 109     Float_t dZEMCorrHG=0.;
 110     Float_t tZN2CorrHG[]={0.,0.,0.,0.,0.}, tZP2CorrHG[]={0.,0.,0.,0.,0.};
 111     Float_t tZN1CorrLG[]={0.,0.,0.,0.,0.}, tZP1CorrLG[]={0.,0.,0.,0.,0.};
 112     Float_t dZEMCorrLG=0.;
 113     Float_t tZN2CorrLG[]={0.,0.,0.,0.,0.}, tZP2CorrLG[]={0.,0.,0.,0.,0.};
 114
 115   for (Int_t iDigit = 0; iDigit < digitsTree->GetEntries(); iDigit++) {
 116     digitsTree->GetEntry(iDigit);
 117     if (!pdigit) continue;
 118
 119     Int_t det = digit.GetSector(0);
 120     Int_t quad = digit.GetSector(1);
 121     Int_t pedindex;
 122     //
 123     if(det == 1){ // *** ZN1
 124        pedindex = quad;
 125        tZN1CorrHG[quad] = (Float_t) (digit.GetADCValue(0)-meanPed[pedindex]);
 126        tZN1CorrLG[quad] = (Float_t) (digit.GetADCValue(1)-meanPed[pedindex+5]);
 127     }
 128     else if(det == 2){ // *** ZP1
 129        pedindex = quad+10;
 130        tZP1CorrHG[quad] = (Float_t) (digit.GetADCValue(0)-meanPed[pedindex]);
 131        tZP1CorrLG[quad] = (Float_t) (digit.GetADCValue(1)-meanPed[pedindex+5]);
 132     }
 133     else if(det == 3){
 134        if(quad == 1){       // *** ZEM1
 135          pedindex = quad+20;
 136          dZEMCorrHG += (Float_t) (digit.GetADCValue(0)-meanPed[pedindex]);
 137          dZEMCorrLG += (Float_t) (digit.GetADCValue(1)-meanPed[pedindex+2]);
 138        }
 139        else if(quad == 2){  // *** ZEM1
 140          pedindex = quad+21;
 141          dZEMCorrHG += (Float_t) (digit.GetADCValue(0)-meanPed[pedindex]);
 142          dZEMCorrLG += (Float_t) (digit.GetADCValue(1)-meanPed[pedindex+2]);
 143        }
 144     }
 145     else if(det == 4){  // *** ZN2
 146        pedindex = quad+24;
 147        tZN2CorrHG[quad] = (Float_t) (digit.GetADCValue(0)-meanPed[pedindex]);
 148        tZN2CorrLG[quad] = (Float_t) (digit.GetADCValue(1)-meanPed[pedindex+5]);
 149     }
 150     else if(det == 5){  // *** ZP2
 151        pedindex = quad+34;
 152        tZP2CorrHG[quad] = (Float_t) (digit.GetADCValue(0)-meanPed[pedindex]);
 153        tZP2CorrLG[quad] = (Float_t) (digit.GetADCValue(1)-meanPed[pedindex+5]);
 154     }
 155   }
 156
 157   // reconstruct the event
 158     ReconstructEvent(clustersTree, tZN1CorrHG, tZP1CorrHG, tZN2CorrHG,
 159         tZP2CorrHG, tZN1CorrLG, tZP1CorrLG, tZN2CorrLG,
 160         tZP2CorrLG, dZEMCorrHG);
 161
 162 }
 163
 164 //_____________________________________________________________________________
 165 void AliZDCReconstructor::Reconstruct(AliRawReader* rawReader, TTree* clustersTree) const
 166 {
 167   // *** ZDC raw data reconstruction
 168   // Works on the current event
 169
 170   // Retrieving calibration data
 171   Float_t meanPed[47];
 172   for(Int_t jj=0; jj<47; jj++) meanPed[jj] = fCalibData->GetMeanPed(jj);
 173
 174   rawReader->Reset();
 175
 176   // loop over raw data rawDatas
 177   Float_t tZN1CorrHG[]={0.,0.,0.,0.,0.}, tZP1CorrHG[]={0.,0.,0.,0.,0.};
 178   Float_t dZEMCorrHG=0.;
 179   Float_t tZN2CorrHG[]={0.,0.,0.,0.,0.}, tZP2CorrHG[]={0.,0.,0.,0.,0.};
 180   Float_t tZN1CorrLG[]={0.,0.,0.,0.,0.}, tZP1CorrLG[]={0.,0.,0.,0.,0.};
 181   Float_t dZEMCorrLG=0.;
 182   Float_t tZN2CorrLG[]={0.,0.,0.,0.,0.}, tZP2CorrLG[]={0.,0.,0.,0.,0.};
 183   //
 184   AliZDCRawStream rawData(rawReader);
 185   while (rawData.Next()) {
 186     if(rawData.IsADCDataWord()){
 187       Int_t det = rawData.GetSector(0);
 188       Int_t quad = rawData.GetSector(1);
 189       Int_t gain = rawData.GetADCGain();
 190       Int_t pedindex;
 191       //
 192       if(det == 1){
 193         pedindex = quad;
 194         if(gain == 0) tZN1CorrHG[quad]  += (Float_t) (rawData.GetADCValue()-meanPed[pedindex]);
 195         else tZN1CorrLG[quad]  += (Float_t) (rawData.GetADCValue()-meanPed[pedindex+5]);
 196       }
 197       else if(det == 2){
 198         pedindex = quad+10;
 199         if(gain == 0) tZP1CorrHG[quad]  += (Float_t) (rawData.GetADCValue()-meanPed[pedindex]);
 200         else tZP1CorrLG[quad]  += (Float_t) (rawData.GetADCValue()-meanPed[pedindex+5]);
 201       }
 202       else if(det == 3){
 203         if(quad==1){
 204           pedindex = quad+20;
 205           if(gain == 0) dZEMCorrHG += (Float_t) (rawData.GetADCValue()-meanPed[pedindex]);
 206           else dZEMCorrLG += (Float_t) (rawData.GetADCValue()-meanPed[pedindex+2]);
 207         }
 208         else if(quad==2){
 209           pedindex = rawData.GetSector(1)+21;
 210           if(gain == 0) dZEMCorrHG += (Float_t) (rawData.GetADCValue()-meanPed[pedindex]);
 211           else dZEMCorrLG += (Float_t) (rawData.GetADCValue()-meanPed[pedindex+2]);
 212         }
 213       }
 214       else if(det == 4){
 215         pedindex = rawData.GetSector(1)+24;
 216         if(gain == 0) tZN2CorrHG[quad]  += (Float_t) (rawData.GetADCValue()-meanPed[pedindex]);
 217         else tZN2CorrLG[quad]  += (Float_t) (rawData.GetADCValue()-meanPed[pedindex+2]);
 218       }
 219       else if(det == 5){
 220         pedindex = rawData.GetSector(1)+34;
 221         if(gain == 0) tZP2CorrHG[quad]  += (Float_t) (rawData.GetADCValue()-meanPed[pedindex]);
 222         else tZP2CorrLG[quad]  += (Float_t) (rawData.GetADCValue()-meanPed[pedindex+5]);
 223       }
 224     }
 225   }
 226
 227   // reconstruct the event
 228     ReconstructEvent(clustersTree, tZN1CorrHG, tZP1CorrHG, tZN2CorrHG,
 229         tZP2CorrHG, tZN1CorrLG, tZP1CorrLG, tZN2CorrLG,
 230         tZP2CorrLG, dZEMCorrHG);
 231
 232 }
 233
 234 //_____________________________________________________________________________
 235 void AliZDCReconstructor::ReconstructEvent(TTree *clustersTree,
 236                 Float_t* ZN1ADCCorrHG, Float_t* ZP1ADCCorrHG,
 237                 Float_t* ZN2ADCCorrHG, Float_t* ZP2ADCCorrHG,
 238                 Float_t* ZN1ADCCorrLG, Float_t* ZP1ADCCorrLG,
 239                 Float_t* ZN2ADCCorrLG, Float_t* ZP2ADCCorrLG,
 240                 Float_t corrADCZEMHG) const
 241 {
 242   // ***** Reconstruct one event
 243
 244   // *** RECONSTRUCTION FROM SIMULATED DATA
 245   // It passes trhough the no. of phe which is known from simulations
 246   //  ---      ADCchannel -> photoelectrons
 247   // NB-> PM gain = 10^(5), ADC resolution = 6.4*10^(-7)
 248   // Move to V965 (E.S.,15/09/04) NB-> PM gain = 10^(5), ADC resolution = 8*10^(-7)
 249   //Float_t zn1phe, zp1phe, zemphe, zn2phe, zp2phe, convFactor = 0.08;
 250   //zn1phe  = ZN1Corr/convFactor;
 251   //zp1phe  = ZP1Corr/convFactor;
 252   //zemphe = ZEMCorr/convFactor;
 253   //zn2phe  = ZN2Corr/convFactor;
 254   //zp2phe  = ZP2Corr/convFactor;
 255   ////if AliDebug(1,Form("\n    znphe = %f, zpphe = %f, zemphe = %f\n",znphe, zpphe, zemphe);
 256   //
 257   ////  ---      Energy calibration
 258   //// Conversion factors for hadronic ZDCs goes from phe yield to TRUE
 259   //// incident energy (conversion from GeV to TeV is included); while for EM
 260   //// calos conversion is from light yield to detected energy calculated by
 261   //// GEANT NB -> ZN and ZP conversion factors are constant since incident
 262   //// spectators have all the same energy, ZEM energy is obtained through a
 263   //// fit over the whole range of incident particle energies
 264   //// (obtained with full HIJING simulations)
 265   //Float_t zn1energy, zp1energy, zemenergy, zdc1energy, zn2energy, zp2energy, zdc2energy;
 266   //Float_t zn1phexTeV=329., zp1phexTeV=369., zn2phexTeV=329., zp2phexTeV=369.;
 267   //zn1energy  = zn1phe/zn1phexTeV;
 268   //zp1energy  = zp1phe/zp1phexTeV;
 269   //zdc1energy = zn1energy+zp1energy;
 270   //zn2energy  = zn2phe/zn2phexTeV;
 271   //zp2energy  = zp2phe/zp2phexTeV;
 272   //zdc2energy = zn2energy+zp2energy;
 273   //zemenergy = -4.81+0.3238*zemphe;
 274   //if(zemenergy<0) zemenergy=0;
 275   ////  if AliDebug(1,Form("    znenergy = %f TeV, zpenergy = %f TeV, zdcenergy = %f GeV, "
 276   ////                     "\n          zemenergy = %f TeV\n", znenergy, zpenergy,
 277   ////                     zdcenergy, zemenergy);
 278   ////  if(zdcenergy==0)
 279   ////    if AliDebug(1,Form("\n\n      ###     ATTENZIONE!!! -> ev# %d: znenergy = %f TeV, zpenergy = %f TeV, zdcenergy = %f GeV, "
 280   ////                       " zemenergy = %f TeV\n\n", fMerger->EvNum(), znenergy, zpenergy, zdcenergy, zemenergy);
 281
 282   //
 283   // *** RECONSTRUCTION FROM "REAL" DATA
 284   //
 285   // Retrieving calibration data
 286   Float_t equalCoeffZN1[5], equalCoeffZP1[5], equalCoeffZN2[5], equalCoeffZP2[5];
 287   for(Int_t ji=0; ji<5; ji++){
 288      equalCoeffZN1[ji] = fCalibData->GetZN1EqualCoeff(ji);
 289      equalCoeffZP1[ji] = fCalibData->GetZP1EqualCoeff(ji);
 290      equalCoeffZN2[ji] = fCalibData->GetZN2EqualCoeff(ji);
 291      equalCoeffZP2[ji] = fCalibData->GetZP2EqualCoeff(ji);
 292   }
 293   //
 294   Float_t calibEne[4];
 295   for(Int_t ij=0; ij<4; ij++) calibEne[ij] = fCalibData->GetEnCalib(ij);
 296   //
 297   Float_t endPointZEM = fCalibData->GetZEMEndValue();
 298   Float_t cutFractionZEM = fCalibData->GetZEMCutFraction();
 299   Float_t dZEMSup = fCalibData->GetDZEMSup();
 300   Float_t dZEMInf = fCalibData->GetDZEMInf();
 301   //
 302   Float_t cutValueZEM = endPointZEM*cutFractionZEM;
 303   Float_t supValueZEM = cutValueZEM+(endPointZEM*dZEMSup);
 304   Float_t infValueZEM = cutValueZEM-(endPointZEM*dZEMInf);
 305   //
 306   Float_t maxValEZN1 = fCalibData->GetEZN1MaxValue();
 307   Float_t maxValEZP1 = fCalibData->GetEZP1MaxValue();
 308   Float_t maxValEZDC1 = fCalibData->GetEZDC1MaxValue();
 309   Float_t maxValEZN2 = fCalibData->GetEZN1MaxValue();
 310   Float_t maxValEZP2 = fCalibData->GetEZP1MaxValue();
 311   Float_t maxValEZDC2 = fCalibData->GetEZDC1MaxValue();
 312
 313   // Equalization of detector responses
 314   Float_t equalTowZN1HG[5], equalTowZN2HG[5], equalTowZP1HG[5], equalTowZP2HG[5];
 315   Float_t equalTowZN1LG[5], equalTowZN2LG[5], equalTowZP1LG[5], equalTowZP2LG[5];
 316   for(Int_t gi=0; gi<5; gi++){
 317      equalTowZN1HG[gi] = ZN1ADCCorrHG[gi]*equalCoeffZN1[gi];
 318      equalTowZP1HG[gi] = ZP1ADCCorrHG[gi]*equalCoeffZP1[gi];
 319      equalTowZN2HG[gi] = ZN2ADCCorrHG[gi]*equalCoeffZN2[gi];
 320      equalTowZP2HG[gi] = ZP2ADCCorrHG[gi]*equalCoeffZP2[gi];
 321      //
 322      equalTowZN1LG[gi] = ZN1ADCCorrLG[gi]*equalCoeffZN1[gi];
 323      equalTowZP1LG[gi] = ZP1ADCCorrLG[gi]*equalCoeffZP1[gi];
 324      equalTowZN2LG[gi] = ZN2ADCCorrLG[gi]*equalCoeffZN2[gi];
 325      equalTowZP2LG[gi] = ZP2ADCCorrLG[gi]*equalCoeffZP2[gi];
 326   }
 327
 328   // Energy calibration of detector responses
 329   Float_t calibTowZN1HG[5], calibTowZN2HG[5], calibTowZP1HG[5], calibTowZP2HG[5];
 330   Float_t calibSumZN1HG=0., calibSumZN2HG=0., calibSumZP1HG=0., calibSumZP2HG=0.;
 331   Float_t calibTowZN1LG[5], calibTowZN2LG[5], calibTowZP1LG[5], calibTowZP2LG[5];
 332   Float_t calibSumZN1LG=0., calibSumZN2LG=0., calibSumZ12LG=0., calibSumZP2LG=0.;
 333   for(Int_t gi=0; gi<5; gi++){
 334      calibTowZN1HG[gi] = equalTowZN1HG[gi]*calibEne[0];
 335      calibTowZP1HG[gi] = equalTowZP1HG[gi]*calibEne[1];
 336      calibTowZN2HG[gi] = equalTowZN2HG[gi]*calibEne[2];
 337      calibTowZP2HG[gi] = equalTowZP2HG[gi]*calibEne[3];
 338      calibSumZN1HG += calibTowZN1HG[gi];
 339      calibSumZP1HG += calibTowZP1HG[gi];
 340      calibSumZN2HG += calibTowZN2HG[gi];
 341      calibSumZP2HG += calibTowZP2HG[gi];
 342      //
 343      calibTowZN1LG[gi] = equalTowZN1LG[gi]*calibEne[0];
 344      calibTowZP1LG[gi] = equalTowZP1LG[gi]*calibEne[1];
 345      calibTowZN2LG[gi] = equalTowZN2LG[gi]*calibEne[2];
 346      calibTowZP2LG[gi] = equalTowZP2LG[gi]*calibEne[3];
 347      calibSumZN1LG += calibTowZN1LG[gi];
 348      calibSumZ12LG += calibTowZP1LG[gi];
 349      calibSumZN2LG += calibTowZN2LG[gi];
 350      calibSumZP2LG += calibTowZP2LG[gi];
 351   }
 352
 353   //  ---      Number of detected spectator nucleons
 354   //  *** N.B. -> It works only in Pb-Pb
 355   Int_t nDetSpecNLeft, nDetSpecPLeft, nDetSpecNRight, nDetSpecPRight;
 356   nDetSpecNLeft = (Int_t) (calibSumZN1HG/2.760);
 357   nDetSpecPLeft = (Int_t) (calibSumZP1HG/2.760);
 358   nDetSpecNRight = (Int_t) (calibSumZN2HG/2.760);
 359   nDetSpecPRight = (Int_t) (calibSumZP2HG/2.760);
 360
 361   //  ---      Number of generated spectator nucleons (from HIJING parameterization)
 362   Int_t nGenSpecNLeft=0, nGenSpecPLeft=0, nGenSpecLeft=0;
 363   Int_t nGenSpecNRight=0, nGenSpecPRight=0, nGenSpecRight=0;
 364   Double_t impPar=0.;
 365   //
 366   // *** RECONSTRUCTION FROM SIMULATED DATA
 367   // Cut value for Ezem (GeV)
 368   // ### Results from production  -> 0<b<18 fm (Apr 2002)
 369   /*Float_t eZEMCut = 420.;
 370   Float_t deltaEZEMSup = 690.;
 371   Float_t deltaEZEMInf = 270.;
 372   if(zemenergy > (eZEMCut+deltaEZEMSup)){
 373     nGenSpecNLeft  = (Int_t) (fZNCen->Eval(ZN1CalibSum));
 374     nGenSpecPLeft  = (Int_t) (fZPCen->Eval(ZP1CalibSum));
 375     nGenSpecLeft   = (Int_t) (fZDCCen->Eval(ZN1CalibSum+ZP1CalibSum));
 376     nGenSpecNRight = (Int_t) (fZNCen->Eval(ZN2CalibSum));
 377     nGenSpecPRight = (Int_t) (fZNCen->Eval(ZP2CalibSum));
 378     nGenSpecRight  = (Int_t) (fZNCen->Eval(ZN2CalibSum+ZP2CalibSum));
 379     impPar  = fbCen->Eval(ZN1CalibSum+ZP1CalibSum);
 380   }
 381   else if(zemenergy < (eZEMCut-deltaEZEMInf)){
 382     nGenSpecNLeft = (Int_t) (fZNPer->Eval(ZN1CalibSum));
 383     nGenSpecPLeft = (Int_t) (fZPPer->Eval(ZP1CalibSum));
 384     nGenSpecLeft  = (Int_t) (fZDCPer->Eval(ZN1CalibSum+ZP1CalibSum));
 385     impPar   = fbPer->Eval(ZN1CalibSum+ZP1CalibSum);
 386   }
 387   else if(zemenergy >= (eZEMCut-deltaEZEMInf) && zemenergy <= (eZEMCut+deltaEZEMSup)){
 388     nGenSpecNLeft = (Int_t) (fZEMn->Eval(zemenergy));
 389     nGenSpecPLeft = (Int_t) (fZEMp->Eval(zemenergy));
 390     nGenSpecLeft  = (Int_t)(fZEMsp->Eval(zemenergy));
 391     impPar   =  fZEMb->Eval(zemenergy);
 392   }
 393   // ### Results from production  -> 0<b<18 fm (Apr 2002)
 394   if(ZN1CalibSum>162.)  nGenSpecNLeft = (Int_t) (fZEMn->Eval(zemenergy));
 395   if(ZP1CalibSum>59.75)  nGenSpecPLeft = (Int_t) (fZEMp->Eval(zemenergy));
 396   if(ZN1CalibSum+ZP1CalibSum>221.5) nGenSpecLeft  = (Int_t)(fZEMsp->Eval(zemenergy));
 397   if(ZN1CalibSum+ZP1CalibSum>220.)  impPar    =  fZEMb->Eval(zemenergy);
 398   */
 399   //
 400   //
 401   // *** RECONSTRUCTION FROM REAL DATA
 402   //
 403   if(corrADCZEMHG > supValueZEM){
 404     nGenSpecNLeft  = (Int_t) (fZNCen->Eval(calibSumZN1HG));
 405     nGenSpecPLeft  = (Int_t) (fZPCen->Eval(calibSumZP1HG));
 406     nGenSpecLeft   = (Int_t) (fZDCCen->Eval(calibSumZN1HG+calibSumZP1HG));
 407     nGenSpecNRight = (Int_t) (fZNCen->Eval(calibSumZN2HG));
 408     nGenSpecPRight = (Int_t) (fZNCen->Eval(calibSumZP2HG));
 409     nGenSpecRight  = (Int_t) (fZNCen->Eval(calibSumZN2HG+calibSumZP2HG));
 410     impPar  = fbCen->Eval(calibSumZN1HG+calibSumZP1HG);
 411   }
 412   else if(corrADCZEMHG < infValueZEM){
 413     nGenSpecNLeft = (Int_t) (fZNPer->Eval(calibSumZN1HG));
 414     nGenSpecPLeft = (Int_t) (fZPPer->Eval(calibSumZP1HG));
 415     nGenSpecLeft  = (Int_t) (fZDCPer->Eval(calibSumZN1HG+calibSumZP1HG));
 416     impPar   = fbPer->Eval(calibSumZN1HG+calibSumZP1HG);
 417   }
 418   else if(corrADCZEMHG >= infValueZEM && corrADCZEMHG <= supValueZEM){
 419     nGenSpecNLeft = (Int_t) (fZEMn->Eval(corrADCZEMHG));
 420     nGenSpecPLeft = (Int_t) (fZEMp->Eval(corrADCZEMHG));
 421     nGenSpecLeft  = (Int_t)(fZEMsp->Eval(corrADCZEMHG));
 422     impPar   =  fZEMb->Eval(corrADCZEMHG);
 423   }
 424   //
 425   if(calibSumZN1HG/maxValEZN1>1.)  nGenSpecNLeft = (Int_t) (fZEMn->Eval(corrADCZEMHG));
 426   if(calibSumZP1HG/maxValEZP1>1.)  nGenSpecPLeft = (Int_t) (fZEMp->Eval(corrADCZEMHG));
 427   if((calibSumZN1HG+calibSumZP1HG/maxValEZDC1)>1.){
 428      nGenSpecLeft = (Int_t)(fZEMsp->Eval(corrADCZEMHG));
 429      impPar = fZEMb->Eval(corrADCZEMHG);
 430   }
 431   if(calibSumZN2HG/maxValEZN2>1.)  nGenSpecNRight = (Int_t) (fZEMn->Eval(corrADCZEMHG));
 432   if(calibSumZP2HG/maxValEZP2>1.)  nGenSpecPRight = (Int_t) (fZEMp->Eval(corrADCZEMHG));
 433   if((calibSumZN2HG+calibSumZP2HG/maxValEZDC2)>1.) nGenSpecRight = (Int_t)(fZEMsp->Eval(corrADCZEMHG));
 434   //
 435   if(nGenSpecNLeft>125)    nGenSpecNLeft=125;
 436   else if(nGenSpecNLeft<0) nGenSpecNLeft=0;
 437   if(nGenSpecPLeft>82)     nGenSpecPLeft=82;
 438   else if(nGenSpecPLeft<0) nGenSpecPLeft=0;
 439   if(nGenSpecLeft>207)     nGenSpecLeft=207;
 440   else if(nGenSpecLeft<0)  nGenSpecLeft=0;
 441
 442   //  ---      Number of generated participants (from HIJING parameterization)
 443   Int_t nPart, nPartTotLeft, nPartTotRight;
 444   nPart = 207-nGenSpecNLeft-nGenSpecPLeft;
 445   nPartTotLeft = 207-nGenSpecLeft;
 446   nPartTotRight = 207-nGenSpecRight;
 447
 448   // create the output tree
 449   AliZDCReco reco(calibSumZN1HG, calibSumZP1HG, calibSumZN2HG, calibSumZP2HG,
 450                   calibTowZN1LG, calibTowZN2LG, calibTowZP1LG, calibTowZP2LG,
 451                   corrADCZEMHG,
 452                   nDetSpecNLeft, nDetSpecPLeft, nDetSpecNRight, nDetSpecPRight,
 453                   nGenSpecNLeft, nGenSpecPLeft, nGenSpecLeft, nGenSpecNRight,
 454                   nGenSpecPRight, nGenSpecRight,
 455                   nPartTotLeft, nPartTotRight, impPar);
 456
 457   AliZDCReco* preco = &reco;
 458   const Int_t kBufferSize = 4000;
 459   clustersTree->Branch("ZDC", "AliZDCReco", &preco, kBufferSize);
 460
 461   // write the output tree
 462   clustersTree->Fill();
 463 }
 464
 465 //_____________________________________________________________________________
 466 void AliZDCReconstructor::FillZDCintoESD(TTree *clustersTree, AliESDEvent* esd) const
 467 {
 468   // fill energies and number of participants to the ESD
 469
 470   AliZDCReco reco;
 471   AliZDCReco* preco = &reco;
 472   clustersTree->SetBranchAddress("ZDC", &preco);
 473
 474   clustersTree->GetEntry(0);
 475   esd->SetZDC(reco.GetZN1Energy(), reco.GetZP1Energy(), reco.GetZEMsignal(),
 476               reco.GetZN2Energy(), reco.GetZP2Energy(),
 477               reco.GetNPartLeft());
 478   //
 479   /*Double_t tZN1Ene[4], tZN2Ene[4];
 480   for(Int_t i=0; i<4; i++){
 481      tZN1Ene[i] = reco.GetZN1EnTow(i);
 482      tZN2Ene[i] = reco.GetZN2EnTow(i);
 483   }*/
 484   //esd->SetZN1TowerEnergy(tZN1Ene);
 485   //esd->SetZN2TowerEnergy(tZN2Ene);
 486
 487 }
 488
 489 //_____________________________________________________________________________
 490 AliCDBStorage* AliZDCReconstructor::SetStorage(const char *uri)
 491 {
 492   // Setting the storage
 493
 494   Bool_t deleteManager = kFALSE;
 495
 496   AliCDBManager *manager = AliCDBManager::Instance();
 497   AliCDBStorage *defstorage = manager->GetDefaultStorage();
 498
 499   if(!defstorage || !(defstorage->Contains("ZDC"))){
 500      AliWarning("No default storage set or default storage doesn't contain ZDC!");
 501      manager->SetDefaultStorage(uri);
 502      deleteManager = kTRUE;
 503   }
 504
 505   AliCDBStorage *storage = manager->GetDefaultStorage();
 506
 507   if(deleteManager){
 508     AliCDBManager::Instance()->UnsetDefaultStorage();
 509     defstorage = 0;   // the storage is killed by AliCDBManager::Instance()->Destroy()
 510   }
 511
 512   return storage;
 513 }
 514
 515 //_____________________________________________________________________________
 516 AliZDCCalibData* AliZDCReconstructor::GetCalibData() const
 517 {
 518
 519   // Getting calibration object for ZDC set
 520
 521   AliCDBEntry  *entry = AliCDBManager::Instance()->Get("ZDC/Calib/Data");
 522   if(!entry) AliFatal("No calibration data loaded!");
 523
 524   AliZDCCalibData *calibdata = dynamic_cast<AliZDCCalibData*>  (entry->GetObject());
 525   if(!calibdata)  AliFatal("Wrong calibration object in calibration  file!");
 526
 527   return calibdata;
 528 }