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 ZN1TowCorrHG[5], ZP1TowCorrHG[5], ZEMCorrHG=0.,
 109             ZN2TowCorrHG[5], ZP2TowCorrHG[5];
 110     Float_t ZN1TowCorrLG[5], ZP1TowCorrLG[5], ZEMCorrLG=0.,
 111             ZN2TowCorrLG[5], ZP2TowCorrLG[5];
 112
 113   for (Int_t iDigit = 0; iDigit < digitsTree->GetEntries(); iDigit++) {
 114     digitsTree->GetEntry(iDigit);
 115     if (!pdigit) continue;
 116
 117     Int_t det = digit.GetSector(0);
 118     Int_t quad = digit.GetSector(1);
 119     Int_t pedindex;
 120     //
 121     if(det == 1){ // *** ZN1
 122        pedindex = quad;
 123        ZN1TowCorrHG[quad] = (Float_t) (digit.GetADCValue(0)-meanPed[pedindex]);
 124        ZN1TowCorrLG[quad] = (Float_t) (digit.GetADCValue(1)-meanPed[pedindex+5]);
 125     }
 126     else if(det == 2){ // *** ZP1
 127        pedindex = quad+10;
 128        ZP1TowCorrHG[quad] = (Float_t) (digit.GetADCValue(0)-meanPed[pedindex]);
 129        ZP1TowCorrLG[quad] = (Float_t) (digit.GetADCValue(1)-meanPed[pedindex+5]);
 130     }
 131     else if(det == 3){
 132        if(quad == 1){       // *** ZEM1
 133          pedindex = quad+20;
 134          ZEMCorrHG += (Float_t) (digit.GetADCValue(0)-meanPed[pedindex]);
 135          ZEMCorrLG += (Float_t) (digit.GetADCValue(1)-meanPed[pedindex+2]);
 136        }
 137        else if(quad == 2){  // *** ZEM1
 138          pedindex = quad+21;
 139          ZEMCorrHG += (Float_t) (digit.GetADCValue(0)-meanPed[pedindex]);
 140          ZEMCorrLG += (Float_t) (digit.GetADCValue(1)-meanPed[pedindex+2]);
 141        }
 142     }
 143     else if(det == 4){  // *** ZN2
 144        pedindex = quad+24;
 145        ZN2TowCorrHG[quad] = (Float_t) (digit.GetADCValue(0)-meanPed[pedindex]);
 146        ZN2TowCorrLG[quad] = (Float_t) (digit.GetADCValue(1)-meanPed[pedindex+5]);
 147     }
 148     else if(det == 5){  // *** ZP2
 149        pedindex = quad+34;
 150        ZP2TowCorrHG[quad] = (Float_t) (digit.GetADCValue(0)-meanPed[pedindex]);
 151        ZP2TowCorrLG[quad] = (Float_t) (digit.GetADCValue(1)-meanPed[pedindex+5]);
 152     }
 153   }
 154
 155   // reconstruct the event
 156     ReconstructEvent(clustersTree, ZN1TowCorrHG, ZP1TowCorrHG, ZN2TowCorrHG,
 157         ZP2TowCorrHG, ZN1TowCorrLG, ZP1TowCorrLG, ZN2TowCorrLG,
 158         ZP2TowCorrLG, ZEMCorrHG);
 159
 160 }
 161
 162 //_____________________________________________________________________________
 163 void AliZDCReconstructor::Reconstruct(AliRawReader* rawReader, TTree* clustersTree) const
 164 {
 165   // *** ZDC raw data reconstruction
 166   // Works on the current event
 167
 168   // Retrieving calibration data
 169   Float_t meanPed[47];
 170   for(Int_t jj=0; jj<47; jj++) meanPed[jj] = fCalibData->GetMeanPed(jj);
 171
 172   rawReader->Reset();
 173
 174   // loop over raw data rawDatas
 175   Float_t ZN1TowCorrHG[5], ZP1TowCorrHG[5], ZEMCorrHG=0.,
 176           ZN2TowCorrHG[5], ZP2TowCorrHG[5];
 177   Float_t ZN1TowCorrLG[5], ZP1TowCorrLG[5], ZEMCorrLG=0.,
 178           ZN2TowCorrLG[5], ZP2TowCorrLG[5];
 179   //
 180   AliZDCRawStream rawData(rawReader);
 181   while (rawData.Next()) {
 182     if(rawData.IsADCDataWord()){
 183       Int_t det = rawData.GetSector(0);
 184       Int_t quad = rawData.GetSector(1);
 185       Int_t gain = rawData.GetADCGain();
 186       Int_t pedindex;
 187       //
 188       if(det == 1){
 189         pedindex = quad;
 190         if(gain == 0) ZN1TowCorrHG[quad]  += (Float_t) (rawData.GetADCValue()-meanPed[pedindex]);
 191         else ZN1TowCorrLG[quad]  += (Float_t) (rawData.GetADCValue()-meanPed[pedindex+5]);
 192       }
 193       else if(det == 2){
 194         pedindex = quad+10;
 195         if(gain == 0) ZP1TowCorrHG[quad]  += (Float_t) (rawData.GetADCValue()-meanPed[pedindex]);
 196         else ZP1TowCorrLG[quad]  += (Float_t) (rawData.GetADCValue()-meanPed[pedindex+5]);
 197       }
 198       else if(det == 3){
 199         if(quad==1){
 200           pedindex = quad+20;
 201           if(gain == 0) ZEMCorrHG += (Float_t) (rawData.GetADCValue()-meanPed[pedindex]);
 202           else ZEMCorrLG += (Float_t) (rawData.GetADCValue()-meanPed[pedindex+2]);
 203         }
 204         else if(quad==2){
 205           pedindex = rawData.GetSector(1)+21;
 206           if(gain == 0) ZEMCorrHG += (Float_t) (rawData.GetADCValue()-meanPed[pedindex]);
 207           else ZEMCorrLG += (Float_t) (rawData.GetADCValue()-meanPed[pedindex+2]);
 208         }
 209       }
 210       else if(det == 4){
 211         pedindex = rawData.GetSector(1)+24;
 212         if(gain == 0) ZN2TowCorrHG[quad]  += (Float_t) (rawData.GetADCValue()-meanPed[pedindex]);
 213         else ZN2TowCorrLG[quad]  += (Float_t) (rawData.GetADCValue()-meanPed[pedindex+2]);
 214       }
 215       else if(det == 5){
 216         pedindex = rawData.GetSector(1)+34;
 217         if(gain == 0) ZP2TowCorrHG[quad]  += (Float_t) (rawData.GetADCValue()-meanPed[pedindex]);
 218         else ZP2TowCorrLG[quad]  += (Float_t) (rawData.GetADCValue()-meanPed[pedindex+5]);
 219       }
 220     }
 221   }
 222
 223   // reconstruct the event
 224     ReconstructEvent(clustersTree, ZN1TowCorrHG, ZP1TowCorrHG, ZN2TowCorrHG,
 225         ZP2TowCorrHG, ZN1TowCorrLG, ZP1TowCorrLG, ZN2TowCorrLG,
 226         ZP2TowCorrLG, ZEMCorrHG);
 227
 228 }
 229
 230 //_____________________________________________________________________________
 231 void AliZDCReconstructor::ReconstructEvent(TTree *clustersTree,
 232                 Float_t* ZN1ADCCorrHG, Float_t* ZP1ADCCorrHG,
 233                 Float_t* ZN2ADCCorrHG, Float_t* ZP2ADCCorrHG,
 234                 Float_t* ZN1ADCCorrLG, Float_t* ZP1ADCCorrLG,
 235                 Float_t* ZN2ADCCorrLG, Float_t* ZP2ADCCorrLG,
 236                 Float_t ZEMADCCorrHG) const
 237 {
 238   // ***** Reconstruct one event
 239
 240   // *** RECONSTRUCTION FROM SIMULATED DATA
 241   // It passes trhough the no. of phe which is known from simulations
 242   //  ---      ADCchannel -> photoelectrons
 243   // NB-> PM gain = 10^(5), ADC resolution = 6.4*10^(-7)
 244   // Move to V965 (E.S.,15/09/04) NB-> PM gain = 10^(5), ADC resolution = 8*10^(-7)
 245   //Float_t zn1phe, zp1phe, zemphe, zn2phe, zp2phe, convFactor = 0.08;
 246   //zn1phe  = ZN1Corr/convFactor;
 247   //zp1phe  = ZP1Corr/convFactor;
 248   //zemphe = ZEMCorr/convFactor;
 249   //zn2phe  = ZN2Corr/convFactor;
 250   //zp2phe  = ZP2Corr/convFactor;
 251   ////if AliDebug(1,Form("\n    znphe = %f, zpphe = %f, zemphe = %f\n",znphe, zpphe, zemphe);
 252   //
 253   ////  ---      Energy calibration
 254   //// Conversion factors for hadronic ZDCs goes from phe yield to TRUE
 255   //// incident energy (conversion from GeV to TeV is included); while for EM
 256   //// calos conversion is from light yield to detected energy calculated by
 257   //// GEANT NB -> ZN and ZP conversion factors are constant since incident
 258   //// spectators have all the same energy, ZEM energy is obtained through a
 259   //// fit over the whole range of incident particle energies
 260   //// (obtained with full HIJING simulations)
 261   //Float_t zn1energy, zp1energy, zemenergy, zdc1energy, zn2energy, zp2energy, zdc2energy;
 262   //Float_t zn1phexTeV=329., zp1phexTeV=369., zn2phexTeV=329., zp2phexTeV=369.;
 263   //zn1energy  = zn1phe/zn1phexTeV;
 264   //zp1energy  = zp1phe/zp1phexTeV;
 265   //zdc1energy = zn1energy+zp1energy;
 266   //zn2energy  = zn2phe/zn2phexTeV;
 267   //zp2energy  = zp2phe/zp2phexTeV;
 268   //zdc2energy = zn2energy+zp2energy;
 269   //zemenergy = -4.81+0.3238*zemphe;
 270   //if(zemenergy<0) zemenergy=0;
 271   ////  if AliDebug(1,Form("    znenergy = %f TeV, zpenergy = %f TeV, zdcenergy = %f GeV, "
 272   ////                     "\n          zemenergy = %f TeV\n", znenergy, zpenergy,
 273   ////                     zdcenergy, zemenergy);
 274   ////  if(zdcenergy==0)
 275   ////    if AliDebug(1,Form("\n\n      ###     ATTENZIONE!!! -> ev# %d: znenergy = %f TeV, zpenergy = %f TeV, zdcenergy = %f GeV, "
 276   ////                       " zemenergy = %f TeV\n\n", fMerger->EvNum(), znenergy, zpenergy, zdcenergy, zemenergy);
 277
 278   //
 279   // *** RECONSTRUCTION FROM "REAL" DATA
 280   //
 281   // Retrieving calibration data
 282   Float_t ZN1EqualCoeff[5], ZP1EqualCoeff[5], ZN2EqualCoeff[5], ZP2EqualCoeff[5];
 283   for(Int_t ji=0; ji<5; ji++){
 284      ZN1EqualCoeff[ji] = fCalibData->GetZN1EqualCoeff(ji);
 285      ZP1EqualCoeff[ji] = fCalibData->GetZP1EqualCoeff(ji);
 286      ZN2EqualCoeff[ji] = fCalibData->GetZN2EqualCoeff(ji);
 287      ZP2EqualCoeff[ji] = fCalibData->GetZP2EqualCoeff(ji);
 288   }
 289   //
 290   Float_t CalibEne[4];
 291   for(Int_t ij=0; ij<4; ij++) CalibEne[ij] = fCalibData->GetEnCalib(ij);
 292   //
 293   Float_t ZEMEndPoint = fCalibData->GetZEMEndValue();
 294   Float_t ZEMCutFraction = fCalibData->GetZEMCutFraction();
 295   Float_t DZEMSup = fCalibData->GetDZEMSup();
 296   Float_t DZEMInf = fCalibData->GetDZEMInf();
 297   //
 298   Float_t ZEMCutValue = ZEMEndPoint*ZEMCutFraction;
 299   Float_t ZEMSupValue = ZEMCutValue+(ZEMEndPoint*DZEMSup);
 300   Float_t ZEMInfValue = ZEMCutValue-(ZEMEndPoint*DZEMInf);
 301   //
 302   Float_t EZN1MaxVal = fCalibData->GetEZN1MaxValue();
 303   Float_t EZP1MaxVal = fCalibData->GetEZP1MaxValue();
 304   Float_t EZDC1MaxVal = fCalibData->GetEZDC1MaxValue();
 305   Float_t EZN2MaxVal = fCalibData->GetEZN1MaxValue();
 306   Float_t EZP2MaxVal = fCalibData->GetEZP1MaxValue();
 307   Float_t EZDC2MaxVal = fCalibData->GetEZDC1MaxValue();
 308
 309   // Equalization of detector responses
 310   Float_t ZN1EqualTowHG[5], ZN2EqualTowHG[5], ZP1EqualTowHG[5], ZP2EqualTowHG[5];
 311   Float_t ZN1EqualTowLG[5], ZN2EqualTowLG[5], ZP1EqualTowLG[5], ZP2EqualTowLG[5];
 312   for(Int_t gi=0; gi<5; gi++){
 313      ZN1EqualTowHG[gi] = ZN1ADCCorrHG[gi]*ZN1EqualCoeff[gi];
 314      ZP1EqualTowHG[gi] = ZP1ADCCorrHG[gi]*ZP1EqualCoeff[gi];
 315      ZN2EqualTowHG[gi] = ZN2ADCCorrHG[gi]*ZN2EqualCoeff[gi];
 316      ZP2EqualTowHG[gi] = ZP2ADCCorrHG[gi]*ZP2EqualCoeff[gi];
 317      //
 318      ZN1EqualTowLG[gi] = ZN1ADCCorrLG[gi]*ZN1EqualCoeff[gi];
 319      ZP1EqualTowLG[gi] = ZP1ADCCorrLG[gi]*ZP1EqualCoeff[gi];
 320      ZN2EqualTowLG[gi] = ZN2ADCCorrLG[gi]*ZN2EqualCoeff[gi];
 321      ZP2EqualTowLG[gi] = ZP2ADCCorrLG[gi]*ZP2EqualCoeff[gi];
 322   }
 323
 324   // Energy calibration of detector responses
 325   Float_t ZN1CalibTowHG[5], ZN2CalibTowHG[5], ZP1CalibTowHG[5], ZP2CalibTowHG[5];
 326   Float_t ZN1CalibSumHG=0., ZN2CalibSumHG=0., ZP1CalibSumHG=0., ZP2CalibSumHG=0.;
 327   Float_t ZN1CalibTowLG[5], ZN2CalibTowLG[5], ZP1CalibTowLG[5], ZP2CalibTowLG[5];
 328   Float_t ZN1CalibSumLG=0., ZN2CalibSumLG=0., ZP1CalibSumLG=0., ZP2CalibSumLG=0.;
 329   for(Int_t gi=0; gi<5; gi++){
 330      ZN1CalibTowHG[gi] = ZN1EqualTowHG[gi]*CalibEne[0];
 331      ZP1CalibTowHG[gi] = ZP1EqualTowHG[gi]*CalibEne[1];
 332      ZN2CalibTowHG[gi] = ZN2EqualTowHG[gi]*CalibEne[2];
 333      ZP2CalibTowHG[gi] = ZP2EqualTowHG[gi]*CalibEne[3];
 334      ZN1CalibSumHG += ZN1CalibTowHG[gi];
 335      ZP1CalibSumHG += ZP1CalibTowHG[gi];
 336      ZN2CalibSumHG += ZN2CalibTowHG[gi];
 337      ZP2CalibSumHG += ZP2CalibTowHG[gi];
 338      //
 339      ZN1CalibTowLG[gi] = ZN1EqualTowLG[gi]*CalibEne[0];
 340      ZP1CalibTowLG[gi] = ZP1EqualTowLG[gi]*CalibEne[1];
 341      ZN2CalibTowLG[gi] = ZN2EqualTowLG[gi]*CalibEne[2];
 342      ZP2CalibTowLG[gi] = ZP2EqualTowLG[gi]*CalibEne[3];
 343      ZN1CalibSumLG += ZN1CalibTowLG[gi];
 344      ZP1CalibSumLG += ZP1CalibTowLG[gi];
 345      ZN2CalibSumLG += ZN2CalibTowLG[gi];
 346      ZP2CalibSumLG += ZP2CalibTowLG[gi];
 347   }
 348
 349   //  ---      Number of detected spectator nucleons
 350   //  *** N.B. -> It works only in Pb-Pb
 351   Int_t nDetSpecNLeft, nDetSpecPLeft, nDetSpecNRight, nDetSpecPRight;
 352   nDetSpecNLeft = (Int_t) (ZN1CalibSumHG/2.760);
 353   nDetSpecPLeft = (Int_t) (ZP1CalibSumHG/2.760);
 354   nDetSpecNRight = (Int_t) (ZN2CalibSumHG/2.760);
 355   nDetSpecPRight = (Int_t) (ZP2CalibSumHG/2.760);
 356
 357   //  ---      Number of generated spectator nucleons (from HIJING parameterization)
 358   Int_t nGenSpecNLeft=0, nGenSpecPLeft=0, nGenSpecLeft=0;
 359   Int_t nGenSpecNRight=0, nGenSpecPRight=0, nGenSpecRight=0;
 360   Double_t impPar=0.;
 361   //
 362   // *** RECONSTRUCTION FROM SIMULATED DATA
 363   // Cut value for Ezem (GeV)
 364   // ### Results from production  -> 0<b<18 fm (Apr 2002)
 365   /*Float_t eZEMCut = 420.;
 366   Float_t deltaEZEMSup = 690.;
 367   Float_t deltaEZEMInf = 270.;
 368   if(zemenergy > (eZEMCut+deltaEZEMSup)){
 369     nGenSpecNLeft  = (Int_t) (fZNCen->Eval(ZN1CalibSum));
 370     nGenSpecPLeft  = (Int_t) (fZPCen->Eval(ZP1CalibSum));
 371     nGenSpecLeft   = (Int_t) (fZDCCen->Eval(ZN1CalibSum+ZP1CalibSum));
 372     nGenSpecNRight = (Int_t) (fZNCen->Eval(ZN2CalibSum));
 373     nGenSpecPRight = (Int_t) (fZNCen->Eval(ZP2CalibSum));
 374     nGenSpecRight  = (Int_t) (fZNCen->Eval(ZN2CalibSum+ZP2CalibSum));
 375     impPar  = fbCen->Eval(ZN1CalibSum+ZP1CalibSum);
 376   }
 377   else if(zemenergy < (eZEMCut-deltaEZEMInf)){
 378     nGenSpecNLeft = (Int_t) (fZNPer->Eval(ZN1CalibSum));
 379     nGenSpecPLeft = (Int_t) (fZPPer->Eval(ZP1CalibSum));
 380     nGenSpecLeft  = (Int_t) (fZDCPer->Eval(ZN1CalibSum+ZP1CalibSum));
 381     impPar   = fbPer->Eval(ZN1CalibSum+ZP1CalibSum);
 382   }
 383   else if(zemenergy >= (eZEMCut-deltaEZEMInf) && zemenergy <= (eZEMCut+deltaEZEMSup)){
 384     nGenSpecNLeft = (Int_t) (fZEMn->Eval(zemenergy));
 385     nGenSpecPLeft = (Int_t) (fZEMp->Eval(zemenergy));
 386     nGenSpecLeft  = (Int_t)(fZEMsp->Eval(zemenergy));
 387     impPar   =  fZEMb->Eval(zemenergy);
 388   }
 389   // ### Results from production  -> 0<b<18 fm (Apr 2002)
 390   if(ZN1CalibSum>162.)  nGenSpecNLeft = (Int_t) (fZEMn->Eval(zemenergy));
 391   if(ZP1CalibSum>59.75)  nGenSpecPLeft = (Int_t) (fZEMp->Eval(zemenergy));
 392   if(ZN1CalibSum+ZP1CalibSum>221.5) nGenSpecLeft  = (Int_t)(fZEMsp->Eval(zemenergy));
 393   if(ZN1CalibSum+ZP1CalibSum>220.)  impPar    =  fZEMb->Eval(zemenergy);
 394   */
 395   //
 396   //
 397   // *** RECONSTRUCTION FROM REAL DATA
 398   //
 399   if(ZEMADCCorrHG > ZEMSupValue){
 400     nGenSpecNLeft  = (Int_t) (fZNCen->Eval(ZN1CalibSumHG));
 401     nGenSpecPLeft  = (Int_t) (fZPCen->Eval(ZP1CalibSumHG));
 402     nGenSpecLeft   = (Int_t) (fZDCCen->Eval(ZN1CalibSumHG+ZP1CalibSumHG));
 403     nGenSpecNRight = (Int_t) (fZNCen->Eval(ZN2CalibSumHG));
 404     nGenSpecPRight = (Int_t) (fZNCen->Eval(ZP2CalibSumHG));
 405     nGenSpecRight  = (Int_t) (fZNCen->Eval(ZN2CalibSumHG+ZP2CalibSumHG));
 406     impPar  = fbCen->Eval(ZN1CalibSumHG+ZP1CalibSumHG);
 407   }
 408   else if(ZEMADCCorrHG < ZEMInfValue){
 409     nGenSpecNLeft = (Int_t) (fZNPer->Eval(ZN1CalibSumHG));
 410     nGenSpecPLeft = (Int_t) (fZPPer->Eval(ZP1CalibSumHG));
 411     nGenSpecLeft  = (Int_t) (fZDCPer->Eval(ZN1CalibSumHG+ZP1CalibSumHG));
 412     impPar   = fbPer->Eval(ZN1CalibSumHG+ZP1CalibSumHG);
 413   }
 414   else if(ZEMADCCorrHG >= ZEMInfValue && ZEMADCCorrHG <= ZEMSupValue){
 415     nGenSpecNLeft = (Int_t) (fZEMn->Eval(ZEMADCCorrHG));
 416     nGenSpecPLeft = (Int_t) (fZEMp->Eval(ZEMADCCorrHG));
 417     nGenSpecLeft  = (Int_t)(fZEMsp->Eval(ZEMADCCorrHG));
 418     impPar   =  fZEMb->Eval(ZEMADCCorrHG);
 419   }
 420   //
 421   if(ZN1CalibSumHG/EZN1MaxVal>1.)  nGenSpecNLeft = (Int_t) (fZEMn->Eval(ZEMADCCorrHG));
 422   if(ZP1CalibSumHG/EZP1MaxVal>1.)  nGenSpecPLeft = (Int_t) (fZEMp->Eval(ZEMADCCorrHG));
 423   if((ZN1CalibSumHG+ZP1CalibSumHG/EZDC1MaxVal)>1.){
 424      nGenSpecLeft = (Int_t)(fZEMsp->Eval(ZEMADCCorrHG));
 425      impPar = fZEMb->Eval(ZEMADCCorrHG);
 426   }
 427   if(ZN2CalibSumHG/EZN2MaxVal>1.)  nGenSpecNRight = (Int_t) (fZEMn->Eval(ZEMADCCorrHG));
 428   if(ZP2CalibSumHG/EZP2MaxVal>1.)  nGenSpecPRight = (Int_t) (fZEMp->Eval(ZEMADCCorrHG));
 429   if((ZN2CalibSumHG+ZP2CalibSumHG/EZDC2MaxVal)>1.) nGenSpecRight = (Int_t)(fZEMsp->Eval(ZEMADCCorrHG));
 430   //
 431   if(nGenSpecNLeft>125)    nGenSpecNLeft=125;
 432   else if(nGenSpecNLeft<0) nGenSpecNLeft=0;
 433   if(nGenSpecPLeft>82)     nGenSpecPLeft=82;
 434   else if(nGenSpecPLeft<0) nGenSpecPLeft=0;
 435   if(nGenSpecLeft>207)     nGenSpecLeft=207;
 436   else if(nGenSpecLeft<0)  nGenSpecLeft=0;
 437
 438   //  ---      Number of generated participants (from HIJING parameterization)
 439   Int_t nPart, nPartTotLeft, nPartTotRight;
 440   nPart = 207-nGenSpecNLeft-nGenSpecPLeft;
 441   nPartTotLeft = 207-nGenSpecLeft;
 442   nPartTotRight = 207-nGenSpecRight;
 443
 444   // create the output tree
 445   AliZDCReco reco(ZN1CalibSumHG, ZP1CalibSumHG, ZN2CalibSumHG, ZP2CalibSumHG,
 446                   ZN1CalibTowLG, ZN2CalibTowLG, ZP1CalibTowLG, ZP2CalibTowLG,
 447                   ZEMADCCorrHG,
 448                   nDetSpecNLeft, nDetSpecPLeft, nDetSpecNRight, nDetSpecPRight,
 449                   nGenSpecNLeft, nGenSpecPLeft, nGenSpecLeft, nGenSpecNRight,
 450                   nGenSpecPRight, nGenSpecRight,
 451                   nPartTotLeft, nPartTotRight, impPar);
 452
 453   AliZDCReco* preco = &reco;
 454   const Int_t kBufferSize = 4000;
 455   clustersTree->Branch("ZDC", "AliZDCReco", &preco, kBufferSize);
 456
 457   // write the output tree
 458   clustersTree->Fill();
 459 }
 460
 461 //_____________________________________________________________________________
 462 void AliZDCReconstructor::FillZDCintoESD(TTree *clustersTree, AliESDEvent* esd) const
 463 {
 464   // fill energies and number of participants to the ESD
 465
 466   AliZDCReco reco;
 467   AliZDCReco* preco = &reco;
 468   clustersTree->SetBranchAddress("ZDC", &preco);
 469
 470   clustersTree->GetEntry(0);
 471   esd->SetZDC(reco.GetZN1Energy(), reco.GetZP1Energy(), reco.GetZEMsignal(),
 472               reco.GetZN2Energy(), reco.GetZP2Energy(),
 473               reco.GetNPartLeft());
 474 }
 475
 476 //_____________________________________________________________________________
 477 AliCDBStorage* AliZDCReconstructor::SetStorage(const char *uri)
 478 {
 479   // Setting the storage
 480
 481   Bool_t deleteManager = kFALSE;
 482
 483   AliCDBManager *manager = AliCDBManager::Instance();
 484   AliCDBStorage *defstorage = manager->GetDefaultStorage();
 485
 486   if(!defstorage || !(defstorage->Contains("ZDC"))){
 487      AliWarning("No default storage set or default storage doesn't contain ZDC!");
 488      manager->SetDefaultStorage(uri);
 489      deleteManager = kTRUE;
 490   }
 491
 492   AliCDBStorage *storage = manager->GetDefaultStorage();
 493
 494   if(deleteManager){
 495     AliCDBManager::Instance()->UnsetDefaultStorage();
 496     defstorage = 0;   // the storage is killed by AliCDBManager::Instance()->Destroy()
 497   }
 498
 499   return storage;
 500 }
 501
 502 //_____________________________________________________________________________
 503 AliZDCCalibData* AliZDCReconstructor::GetCalibData() const
 504 {
 505
 506   // Getting calibration object for ZDC set
 507
 508   AliCDBEntry  *entry = AliCDBManager::Instance()->Get("ZDC/Calib/Data");
 509   if(!entry) AliFatal("No calibration data loaded!");
 510
 511   AliZDCCalibData *calibdata = dynamic_cast<AliZDCCalibData*>  (entry->GetObject());
 512   if(!calibdata)  AliFatal("Wrong calibration object in calibration  file!");
 513
 514   return calibdata;
 515 }