1 /**************************************************************************
2 * Copyright(c) 1998-1999, ALICE Experiment at CERN, All rights reserved. *
4 * Author: The ALICE Off-line Project. *
5 * Contributors are mentioned in the code where appropriate. *
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 **************************************************************************/
18 ///////////////////////////////////////////////////////////////////////////////
20 // ************** Class for ZDC reconstruction ************** //
21 // Author: Chiara.Oppedisano@to.infn.it //
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) //
27 ///////////////////////////////////////////////////////////////////////////////
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 #include "AliLHCClockPhase.h"
53 ClassImp(AliZDCReconstructor)
54 AliZDCRecoParam *AliZDCReconstructor::fgRecoParam=0; //reconstruction parameters
55 AliZDCMBCalib *AliZDCReconstructor::fgMBCalibData=0; //calibration parameters for A-A reconstruction
57 //_____________________________________________________________________________
58 AliZDCReconstructor:: AliZDCReconstructor() :
59 fPedData(GetPedestalData()),
60 fEnCalibData(GetEnergyCalibData()),
61 fTowCalibData(GetTowerCalibData()),
65 fIsCalibrationMB(kFALSE),
71 // **** Default constructor
75 //_____________________________________________________________________________
76 AliZDCReconstructor::~AliZDCReconstructor()
79 // if(fgRecoParam) delete fgRecoParam;
80 if(fPedData) delete fPedData;
81 if(fEnCalibData) delete fEnCalibData;
82 if(fTowCalibData) delete fTowCalibData;
83 if(fgMBCalibData) delete fgMBCalibData;
84 if(fESDZDC) delete fESDZDC;
87 //____________________________________________________________________________
88 void AliZDCReconstructor::Init()
90 // Setting reconstruction parameters
92 TString runType = GetRunInfo()->GetRunType();
93 if((runType.CompareTo("CALIBRATION_MB")) == 0){
94 fIsCalibrationMB = kTRUE;
97 TString beamType = GetRunInfo()->GetBeamType();
98 // This is a temporary solution to allow reconstruction in tests without beam
99 if(((beamType.CompareTo("UNKNOWN"))==0) &&
100 ((runType.CompareTo("PHYSICS"))==0 || (runType.CompareTo("CALIBRATION_BC"))==0)){
103 /*else if((beamType.CompareTo("UNKNOWN"))==0){
104 AliError("\t UNKNOWN beam type\n");
108 fBeamEnergy = GetRunInfo()->GetBeamEnergy();
109 if(fBeamEnergy<0.01) AliWarning(" Beam energy value missing -> E_beam = 0");
111 if(((beamType.CompareTo("pp"))==0) || ((beamType.CompareTo("p-p"))==0)
112 ||((beamType.CompareTo("PP"))==0) || ((beamType.CompareTo("P-P"))==0)){
115 else if((beamType.CompareTo("A-A")) == 0 || (beamType.CompareTo("AA")) == 0){
117 if(!fgRecoParam) fgRecoParam = const_cast<AliZDCRecoParam*>(GetRecoParam());
119 fgRecoParam->SetGlauberMCDist(fBeamEnergy);
123 AliCDBEntry *entry = AliCDBManager::Instance()->Get("GRP/Calib/LHCClockPhase");
124 if (!entry) AliFatal("LHC clock-phase shift is not found in OCDB !");
125 AliLHCClockPhase *phaseLHC = (AliLHCClockPhase*)entry->GetObject();
126 // 4/2/2011 According to A. Di Mauro BEAM1 measurement is more reliable
127 // than BEAM2 and therefore also than the average of the 2
128 fMeanPhase = phaseLHC->GetMeanPhaseB1();
130 if(fIsCalibrationMB==kFALSE)
131 printf("\n\n ***** ZDC reconstruction initialized for %s @ %1.0f + %1.0f GeV *****\n\n",
132 beamType.Data(), fBeamEnergy, fBeamEnergy);
134 // if EMD calibration run NO ENERGY CALIBRATION should be performed
135 // pp-like reconstruction must be performed (E cailb. coeff. = 1)
136 if((runType.CompareTo("CALIBRATION_EMD")) == 0){
141 fESDZDC = new AliESDZDC();
146 //____________________________________________________________________________
147 void AliZDCReconstructor::Init(TString beamType, Float_t beamEnergy)
149 // Setting reconstruction mode
150 // Needed to work in the HLT framework
152 fIsCalibrationMB = kFALSE;
154 fBeamEnergy = beamEnergy;
156 if(((beamType.CompareTo("pp"))==0) || ((beamType.CompareTo("p-p"))==0)
157 ||((beamType.CompareTo("PP"))==0) || ((beamType.CompareTo("P-P"))==0)){
160 else if((beamType.CompareTo("A-A")) == 0 || (beamType.CompareTo("AA")) == 0){
162 if(!fgRecoParam) fgRecoParam = const_cast<AliZDCRecoParam*>(GetRecoParam());
163 if( fgRecoParam ) fgRecoParam->SetGlauberMCDist(fBeamEnergy);
166 AliCDBEntry *entry = AliCDBManager::Instance()->Get("GRP/Calib/LHCClockPhase");
167 if (!entry) AliFatal("LHC clock-phase shift is not found in OCDB !");
168 AliLHCClockPhase *phaseLHC = (AliLHCClockPhase*)entry->GetObject();
169 fMeanPhase = phaseLHC->GetMeanPhase();
171 fESDZDC = new AliESDZDC();
173 printf("\n\n ***** ZDC reconstruction initialized for %s @ %1.0f + %1.0f GeV *****\n\n",
174 beamType.Data(), fBeamEnergy, fBeamEnergy);
178 //_____________________________________________________________________________
179 void AliZDCReconstructor::Reconstruct(TTree* digitsTree, TTree* clustersTree) const
181 // *** Local ZDC reconstruction for digits
182 // Works on the current event
184 // Retrieving calibration data
185 // Parameters for mean value pedestal subtraction
187 Float_t meanPed[2*kNch];
188 for(Int_t jj=0; jj<2*kNch; jj++) meanPed[jj] = fPedData->GetMeanPed(jj);
189 // Parameters pedestal subtraction through correlation with out-of-time signals
190 Float_t corrCoeff0[2*kNch], corrCoeff1[2*kNch];
191 for(Int_t jj=0; jj<2*kNch; jj++){
192 corrCoeff0[jj] = fPedData->GetPedCorrCoeff0(jj);
193 corrCoeff1[jj] = fPedData->GetPedCorrCoeff1(jj);
198 AliZDCDigit* pdigit = &digit;
199 digitsTree->SetBranchAddress("ZDC", &pdigit);
200 //printf("\n\t # of digits in tree: %d\n",(Int_t) digitsTree->GetEntries());
203 Float_t tZN1Corr[10], tZP1Corr[10], tZN2Corr[10], tZP2Corr[10];
204 Float_t dZEM1Corr[2], dZEM2Corr[2], sPMRef1[2], sPMRef2[2];
205 for(Int_t i=0; i<10; i++){
206 tZN1Corr[i] = tZP1Corr[i] = tZN2Corr[i] = tZP2Corr[i] = 0.;
207 if(i<2) dZEM1Corr[i] = dZEM2Corr[i] = sPMRef1[i] = sPMRef2[i] = 0.;
210 Int_t digNentries = digitsTree->GetEntries();
211 Float_t ootDigi[kNch]; Int_t i=0;
212 // -- Reading out-of-time signals (last kNch entries) for current event
214 for(Int_t iDigit=kNch; iDigit<digNentries; iDigit++){
215 if(i<=kNch) ootDigi[i] = digitsTree->GetEntry(iDigit);
216 else AliWarning(" Can't read more out of time values: index>kNch !!!\n");
221 for(Int_t iDigit=0; iDigit<(digNentries/2); iDigit++) {
222 digitsTree->GetEntry(iDigit);
223 if (!pdigit) continue;
225 Int_t det = digit.GetSector(0);
226 Int_t quad = digit.GetSector(1);
228 Float_t ped2SubHg=0., ped2SubLg=0.;
230 if(det==1) pedindex = quad;
231 else if(det==2) pedindex = quad+5;
232 else if(det==3) pedindex = quad+9;
233 else if(det==4) pedindex = quad+12;
234 else if(det==5) pedindex = quad+17;
236 else pedindex = (det-1)/3+22;
239 ped2SubHg = meanPed[pedindex];
240 ped2SubLg = meanPed[pedindex+kNch];
242 else if(fPedSubMode==1){
243 ped2SubHg = corrCoeff1[pedindex]*ootDigi[pedindex]+corrCoeff0[pedindex];
244 ped2SubLg = corrCoeff1[pedindex+kNch]*ootDigi[pedindex+kNch]+corrCoeff0[pedindex+kNch];
247 if(quad != 5){ // ZDC (not reference PTMs!)
248 if(det == 1){ // *** ZNC
249 tZN1Corr[quad] = (Float_t) (digit.GetADCValue(0)-ped2SubHg);
250 tZN1Corr[quad+5] = (Float_t) (digit.GetADCValue(1)-ped2SubLg);
251 if(tZN1Corr[quad]<0.) tZN1Corr[quad] = 0.;
252 if(tZN1Corr[quad+5]<0.) tZN1Corr[quad+5] = 0.;
254 else if(det == 2){ // *** ZP1
255 tZP1Corr[quad] = (Float_t) (digit.GetADCValue(0)-ped2SubHg);
256 tZP1Corr[quad+5] = (Float_t) (digit.GetADCValue(1)-ped2SubLg);
257 if(tZP1Corr[quad]<0.) tZP1Corr[quad] = 0.;
258 if(tZP1Corr[quad+5]<0.) tZP1Corr[quad+5] = 0.;
261 if(quad == 1){ // *** ZEM1
262 dZEM1Corr[0] += (Float_t) (digit.GetADCValue(0)-ped2SubHg);
263 dZEM1Corr[1] += (Float_t) (digit.GetADCValue(1)-ped2SubLg);
264 if(dZEM1Corr[0]<0.) dZEM1Corr[0] = 0.;
265 if(dZEM1Corr[1]<0.) dZEM1Corr[1] = 0.;
267 else if(quad == 2){ // *** ZEM2
268 dZEM2Corr[0] += (Float_t) (digit.GetADCValue(0)-ped2SubHg);
269 dZEM2Corr[1] += (Float_t) (digit.GetADCValue(1)-ped2SubLg);
270 if(dZEM2Corr[0]<0.) dZEM2Corr[0] = 0.;
271 if(dZEM2Corr[1]<0.) dZEM2Corr[1] = 0.;
274 else if(det == 4){ // *** ZN2
275 tZN2Corr[quad] = (Float_t) (digit.GetADCValue(0)-ped2SubHg);
276 tZN2Corr[quad+5] = (Float_t) (digit.GetADCValue(1)-ped2SubLg);
277 if(tZN2Corr[quad]<0.) tZN2Corr[quad] = 0.;
278 if(tZN2Corr[quad+5]<0.) tZN2Corr[quad+5] = 0.;
280 else if(det == 5){ // *** ZP2
281 tZP2Corr[quad] = (Float_t) (digit.GetADCValue(0)-ped2SubHg);
282 tZP2Corr[quad+5] = (Float_t) (digit.GetADCValue(1)-ped2SubLg);
283 if(tZP2Corr[quad]<0.) tZP2Corr[quad] = 0.;
284 if(tZP2Corr[quad+5]<0.) tZP2Corr[quad+5] = 0.;
287 else{ // Reference PMs
289 sPMRef1[0] = (Float_t) (digit.GetADCValue(0)-ped2SubHg);
290 sPMRef1[1] = (Float_t) (digit.GetADCValue(1)-ped2SubLg);
292 if(sPMRef1[0]<0.) sPMRef1[0] = 0.;
293 if(sPMRef2[1]<0.) sPMRef1[1] = 0.;
296 sPMRef2[0] = (Float_t) (digit.GetADCValue(0)-ped2SubHg);
297 sPMRef2[1] = (Float_t) (digit.GetADCValue(1)-ped2SubLg);
299 if(sPMRef2[0]<0.) sPMRef2[0] = 0.;
300 if(sPMRef2[1]<0.) sPMRef2[1] = 0.;
305 /*printf("AliZDCReconstructor: digit #%d det %d quad %d pedHG %1.0f pedLG %1.0f\n",
306 iDigit, det, quad, ped2SubHg, ped2SubLg);
307 printf(" -> pedindex %d\n", pedindex);
308 printf(" HGChain -> RawDig %d DigCorr %1.2f",
309 digit.GetADCValue(0), digit.GetADCValue(0)-ped2SubHg);
310 printf(" LGChain -> RawDig %d DigCorr %1.2f\n",
311 digit.GetADCValue(1), digit.GetADCValue(1)-ped2SubLg);*/
317 for(Int_t jj=0; jj<32; jj++){
319 for(Int_t ii=0; ii<4; ii++) tdc[jj][ii]=0;
322 Int_t evQualityBlock[4] = {1,0,0,0};
323 Int_t triggerBlock[4] = {0,0,0,0};
324 Int_t chBlock[3] = {0,0,0};
327 // reconstruct the event
329 ReconstructEventpp(clustersTree, tZN1Corr, tZP1Corr, tZN2Corr, tZP2Corr,
330 dZEM1Corr, dZEM2Corr, sPMRef1, sPMRef2,
332 evQualityBlock, triggerBlock, chBlock, puBits);
333 else if(fRecoMode==2)
334 ReconstructEventPbPb(clustersTree, tZN1Corr, tZP1Corr, tZN2Corr, tZP2Corr,
335 dZEM1Corr, dZEM2Corr, sPMRef1, sPMRef2,
337 evQualityBlock, triggerBlock, chBlock, puBits);
340 //_____________________________________________________________________________
341 void AliZDCReconstructor::Reconstruct(AliRawReader* rawReader, TTree* clustersTree) const
343 // *** ZDC raw data reconstruction
344 // Works on the current event
346 // Retrieving calibration data
347 // Parameters for pedestal subtraction
349 Float_t meanPed[2*kNch];
350 for(Int_t jj=0; jj<2*kNch; jj++) meanPed[jj] = fPedData->GetMeanPed(jj);
351 // Parameters pedestal subtraction through correlation with out-of-time signals
352 Float_t corrCoeff0[2*kNch], corrCoeff1[2*kNch];
353 for(Int_t jj=0; jj<2*kNch; jj++){
354 corrCoeff0[jj] = fPedData->GetPedCorrCoeff0(jj);
355 corrCoeff1[jj] = fPedData->GetPedCorrCoeff1(jj);
356 //printf(" %d %1.4f %1.4f\n", jj,corrCoeff0[jj],corrCoeff1[jj]);
359 Int_t adcZN1[5], adcZN1oot[5], adcZN1lg[5], adcZN1ootlg[5];
360 Int_t adcZP1[5], adcZP1oot[5], adcZP1lg[5], adcZP1ootlg[5];
361 Int_t adcZN2[5], adcZN2oot[5], adcZN2lg[5], adcZN2ootlg[5];
362 Int_t adcZP2[5], adcZP2oot[5], adcZP2lg[5], adcZP2ootlg[5];
363 Int_t adcZEM[2], adcZEMoot[2], adcZEMlg[2], adcZEMootlg[2];
364 Int_t pmRef[2], pmRefoot[2], pmReflg[2], pmRefootlg[2];
365 for(Int_t ich=0; ich<5; ich++){
366 adcZN1[ich] = adcZN1oot[ich] = adcZN1lg[ich] = adcZN1ootlg[ich] = 0;
367 adcZP1[ich] = adcZP1oot[ich] = adcZP1lg[ich] = adcZP1ootlg[ich] = 0;
368 adcZN2[ich] = adcZN2oot[ich] = adcZN2lg[ich] = adcZN2ootlg[ich] = 0;
369 adcZP2[ich] = adcZP2oot[ich] = adcZP2lg[ich] = adcZP2ootlg[ich] = 0;
371 adcZEM[ich] = adcZEMoot[ich] = adcZEMlg[ich] = adcZEMootlg[ich] = 0;
372 pmRef[ich] = pmRefoot[ich] = pmReflg[ich] = pmRefootlg[ich] = 0;
376 Float_t tZN1Corr[10], tZP1Corr[10], tZN2Corr[10], tZP2Corr[10];
377 Float_t dZEM1Corr[2], dZEM2Corr[2], sPMRef1[2], sPMRef2[2];
378 for(Int_t i=0; i<10; i++){
379 tZN1Corr[i] = tZP1Corr[i] = tZN2Corr[i] = tZP2Corr[i] = 0.;
380 if(i<2) dZEM1Corr[i] = dZEM2Corr[i] = sPMRef1[i] = sPMRef2[i] = 0.;
383 Bool_t isScalerOn=kFALSE;
384 Int_t jsc=0, itdc=0, iprevtdc=-1, ihittdc=0;
385 UInt_t scalerData[32];
386 Int_t tdcData[32][4];
387 for(Int_t k=0; k<32; k++){
389 for(Int_t i=0; i<4; i++) tdcData[k][i]=0;
393 Int_t evQualityBlock[4] = {1,0,0,0};
394 Int_t triggerBlock[4] = {0,0,0,0};
395 Int_t chBlock[3] = {0,0,0};
398 Int_t kFirstADCGeo=0, kLastADCGeo=3, kScalerGeo=8, kZDCTDCGeo=4, kPUGeo=29;
399 //Int_t kTrigScales=30, kTrigHistory=31;
401 // loop over raw data
402 //rawReader->Reset();
403 AliZDCRawStream rawData(rawReader);
404 while(rawData.Next()){
406 // ***************************** Reading ADCs
407 if((rawData.GetADCModule()>=kFirstADCGeo) && (rawData.GetADCModule()<=kLastADCGeo)){
408 //printf(" **** Reading ADC raw data from module %d **** \n",rawData.GetADCModule());
410 if((rawData.IsADCDataWord()) && (rawData.GetNChannelsOn()<48)) chBlock[0] = kTRUE;
411 if((rawData.IsADCDataWord()) && (rawData.IsOverflow() == kTRUE)) chBlock[1] = kTRUE;
412 if((rawData.IsADCDataWord()) && (rawData.IsUnderflow() == kTRUE)) chBlock[2] = kTRUE;
413 if((rawData.IsADCDataWord()) && (rawData.IsADCEventGood() == kTRUE)) evQualityBlock[0] = kTRUE;
415 if((rawData.IsADCDataWord()) && (rawData.IsUnderflow()==kFALSE)
416 && (rawData.IsOverflow()==kFALSE) && (rawData.IsADCEventGood()==kTRUE)){
418 Int_t adcMod = rawData.GetADCModule();
419 Int_t det = rawData.GetSector(0);
420 Int_t quad = rawData.GetSector(1);
421 Int_t gain = rawData.GetADCGain();
424 // Mean pedestal value subtraction -------------------------------------------------------
425 if(fPedSubMode == 0){
426 // **** Pb-Pb data taking 2010 -> subtracting some ch. from correlation ****
427 // Not interested in o.o.t. signals (ADC modules 2, 3)
428 //if(adcMod == 2 || adcMod == 3) continue;
429 if(((det==1 && quad==0) || (det==3))){
431 if(adcMod==0 || adcMod==1){
432 if(gain==0) adcZN1[quad] = rawData.GetADCValue();
433 else adcZN1lg[quad] = rawData.GetADCValue();
435 else if(adcMod==2 || adcMod==3){
436 if(gain==0) adcZN1oot[quad] = rawData.GetADCValue();
437 else adcZN1ootlg[quad] = rawData.GetADCValue();
441 if(adcMod==0 || adcMod==1){
442 if(gain==0) adcZEM[quad-1] = rawData.GetADCValue();
443 else adcZEMlg[quad-1] = rawData.GetADCValue();
445 else if(adcMod==2 || adcMod==3){
446 if(gain==0) adcZEMoot[quad-1] = rawData.GetADCValue();
447 else adcZEMootlg[quad-1] = rawData.GetADCValue();
451 // When oot values are read the ADC modules 2, 3 can be skipped!!!
452 if(adcMod == 2 || adcMod == 3) continue;
454 // *************************************************************************
455 if(quad != 5){ // ZDCs (not reference PTMs)
456 if(det==1 && quad!=0){
458 if(gain == 0) tZN1Corr[quad] += (Float_t) (rawData.GetADCValue()-meanPed[pedindex]);
459 else tZN1Corr[quad+5] += (Float_t) (rawData.GetADCValue()-meanPed[pedindex+kNch]);
463 if(gain == 0) tZP1Corr[quad] += (Float_t) (rawData.GetADCValue()-meanPed[pedindex]);
464 else tZP1Corr[quad+5] += (Float_t) (rawData.GetADCValue()-meanPed[pedindex+kNch]);
469 if(gain == 0) dZEM1Corr[0] += (Float_t) (rawData.GetADCValue()-meanPed[pedindex]);
470 else dZEM1Corr[1] += (Float_t) (rawData.GetADCValue()-meanPed[pedindex+kNch]);
473 if(gain == 0) dZEM2Corr[0] += (Float_t) (rawData.GetADCValue()-meanPed[pedindex]);
474 else dZEM2Corr[1] += (Float_t) (rawData.GetADCValue()-meanPed[pedindex+kNch]);
479 if(gain == 0) tZN2Corr[quad] += (Float_t) (rawData.GetADCValue()-meanPed[pedindex]);
480 else tZN2Corr[quad+5] += (Float_t) (rawData.GetADCValue()-meanPed[pedindex+kNch]);
484 if(gain == 0) tZP2Corr[quad] += (Float_t) (rawData.GetADCValue()-meanPed[pedindex]);
485 else tZP2Corr[quad+5] += (Float_t) (rawData.GetADCValue()-meanPed[pedindex+kNch]);
488 else{ // reference PM
489 pedindex = (det-1)/3 + 22;
491 if(gain==0) sPMRef1[0] += (Float_t) (rawData.GetADCValue()-meanPed[pedindex]);
492 else sPMRef1[1] += (Float_t) (rawData.GetADCValue()-meanPed[pedindex+kNch]);
495 if(gain==0) sPMRef2[0] += (Float_t) (rawData.GetADCValue()-meanPed[pedindex]);
496 else sPMRef2[1] += (Float_t) (rawData.GetADCValue()-meanPed[pedindex+kNch]);
501 printf(" AliZDCReconstructor: det %d quad %d res %d -> Pedestal[%d] %1.0f",
502 det,quad,gain, pedindex, meanPed[pedindex]);
503 printf(" RawADC %d ADCCorr %1.0f\n",
504 rawData.GetADCValue(), rawData.GetADCValue()-meanPed[pedindex]);
506 }// mean pedestal subtraction
507 // Pedestal subtraction from correlation ------------------------------------------------
508 else if(fPedSubMode == 1){
510 if(adcMod==0 || adcMod==1){
511 if(quad != 5){ // signals from ZDCs
513 if(gain==0) adcZN1[quad] = rawData.GetADCValue();
514 else adcZN1lg[quad] = rawData.GetADCValue();
517 if(gain==0) adcZP1[quad] = rawData.GetADCValue();
518 else adcZP1lg[quad] = rawData.GetADCValue();
521 if(gain==0) adcZEM[quad-1] = rawData.GetADCValue();
522 else adcZEMlg[quad-1] = rawData.GetADCValue();
525 if(gain==0) adcZN2[quad] = rawData.GetADCValue();
526 else adcZN2lg[quad] = rawData.GetADCValue();
529 if(gain==0) adcZP2[quad] = rawData.GetADCValue();
530 else adcZP2lg[quad] = rawData.GetADCValue();
533 else{ // signals from reference PM
534 if(gain==0) pmRef[quad-1] = rawData.GetADCValue();
535 else pmReflg[quad-1] = rawData.GetADCValue();
538 // Out-of-time pedestals
539 else if(adcMod==2 || adcMod==3){
540 if(quad != 5){ // signals from ZDCs
542 if(gain==0) adcZN1oot[quad] = rawData.GetADCValue();
543 else adcZN1ootlg[quad] = rawData.GetADCValue();
546 if(gain==0) adcZP1oot[quad] = rawData.GetADCValue();
547 else adcZP1ootlg[quad] = rawData.GetADCValue();
550 if(gain==0) adcZEMoot[quad-1] = rawData.GetADCValue();
551 else adcZEMootlg[quad-1] = rawData.GetADCValue();
554 if(gain==0) adcZN2oot[quad] = rawData.GetADCValue();
555 else adcZN2ootlg[quad] = rawData.GetADCValue();
558 if(gain==0) adcZP2oot[quad] = rawData.GetADCValue();
559 else adcZP2ootlg[quad] = rawData.GetADCValue();
562 else{ // signals from reference PM
563 if(gain==0) pmRefoot[quad-1] = rawData.GetADCValue();
564 else pmRefootlg[quad-1] = rawData.GetADCValue();
567 } // pedestal subtraction from correlation
569 /*printf("\t AliZDCReconstructor: det %d quad %d res %d -> Ped[%d] = %1.0f\n",
570 det,quad,gain, pedindex, meanPed[pedindex]);*/
573 // ***************************** Reading Scaler
574 else if(rawData.GetADCModule()==kScalerGeo){
575 if(rawData.IsScalerWord()==kTRUE && rawData.IsScEventGood()==kTRUE){
577 scalerData[jsc] = rawData.GetTriggerCount();
579 //printf(" Reconstructed VME Scaler: %d %d ",jsc,scalerData[jsc]);
584 // ***************************** Reading ZDC TDC
585 else if(rawData.GetADCModule()==kZDCTDCGeo && rawData.IsZDCTDCDatum()==kTRUE){
586 itdc = rawData.GetChannel();
587 if(itdc==iprevtdc) ihittdc++;
590 tdcData[itdc][ihittdc] = rawData.GetZDCTDCDatum();
592 //printf(" Reconstructed TDC[%d, %d] %d ",itdc, ihittdc, tdcData[itdc][ihittdc]);
594 // ***************************** Reading PU
595 else if(rawData.GetADCModule()==kPUGeo){
596 puBits = rawData.GetDetectorPattern();
598 // ***************************** Reading trigger history
599 else if(rawData.IstriggerHistoryWord()==kTRUE){
600 triggerBlock[0] = rawData.IsCPTInputEMDTrigger();
601 triggerBlock[1] = rawData.IsCPTInputSemiCentralTrigger();
602 triggerBlock[2] = rawData.IsCPTInputCentralTrigger();
603 triggerBlock[3] = rawData.IsCPTInputMBTrigger();
609 for(Int_t t=0; t<5; t++){
610 tZN1Corr[t] = adcZN1[t] - (corrCoeff1[t]*adcZN1oot[t]+corrCoeff0[t]);
611 tZN1Corr[t+5] = adcZN1lg[t] - (corrCoeff1[t+kNch]*adcZN1ootlg[t]+corrCoeff0[t+kNch]);
613 tZP1Corr[t] = adcZP1[t] - (corrCoeff1[t+5]*adcZP1oot[t]+corrCoeff0[t+5]);
614 tZP1Corr[t+5] = adcZP1lg[t] - (corrCoeff1[t+5+kNch]*adcZP1ootlg[t]+corrCoeff0[t+5+kNch]);
616 tZN2Corr[t] = adcZN2[t] - (corrCoeff1[t+12]*adcZN2oot[t]+corrCoeff0[t+12]);
617 tZN2Corr[t+5] = adcZN2lg[t] - (corrCoeff1[t+12+kNch]*adcZN2ootlg[t]+corrCoeff0[t+12+kNch]);
619 tZP2Corr[t] = adcZP2[t] - (corrCoeff1[t+17]*adcZP2oot[t]+corrCoeff0[t+17]);
620 tZP2Corr[t+5] = adcZP2lg[t] - (corrCoeff1[t+17+kNch]*adcZP2ootlg[t]+corrCoeff0[t+17+kNch]);
622 dZEM1Corr[0] = adcZEM[0] - (corrCoeff1[10]*adcZEMoot[0]+corrCoeff0[10]);
623 dZEM1Corr[1] = adcZEMlg[0] - (corrCoeff1[10+kNch]*adcZEMootlg[0]+corrCoeff0[10+kNch]);
624 dZEM2Corr[0] = adcZEM[1] - (corrCoeff1[11]*adcZEMoot[1]+corrCoeff0[11]);
625 dZEM2Corr[1] = adcZEMlg[1] - (corrCoeff1[11+kNch]*adcZEMootlg[1]+corrCoeff0[11+kNch]);
627 sPMRef1[0] = pmRef[0] - (corrCoeff1[22]*pmRefoot[0]+corrCoeff0[22]);
628 sPMRef1[1] = pmReflg[0] - (corrCoeff1[22+kNch]*pmRefootlg[0]+corrCoeff0[22+kNch]);
629 sPMRef2[0] = pmRef[0] - (corrCoeff1[23]*pmRefoot[1]+corrCoeff0[23]);
630 sPMRef2[1] = pmReflg[0] - (corrCoeff1[23+kNch]*pmRefootlg[1]+corrCoeff0[23+kNch]);
633 // **** Pb-Pb data taking 2010 -> subtracting some ch. from correlation ****
634 tZN1Corr[0] = adcZN1[0] - (corrCoeff1[0]*adcZN1oot[0]+corrCoeff0[0]);
635 tZN1Corr[5] = adcZN1lg[0] - (corrCoeff1[kNch]*adcZN1ootlg[0]+corrCoeff0[kNch]);
637 //printf(" adcZN1 %d adcZN1oot %d tZN1Corr %1.2f \n", adcZN1[0],adcZN1oot[0],tZN1Corr[0]);
638 //printf(" adcZN1lg %d adcZN1ootlg %d tZN1Corrlg %1.2f \n", adcZN1lg[0],adcZN1ootlg[0],tZN1Corr[5]);
640 //tZP1Corr[2] = adcZP1[2] - (corrCoeff1[2+5]*adcZP1oot[2]+corrCoeff0[2+5]);
641 //tZP1Corr[2+5] = adcZP1lg[2] - (corrCoeff1[2+5+kNch]*adcZP1ootlg[2]+corrCoeff0[2+5+kNch]);
643 dZEM1Corr[0] = adcZEM[0] - (corrCoeff1[10]*adcZEMoot[0]+corrCoeff0[10]);
644 dZEM1Corr[1] = adcZEMlg[0] - (corrCoeff1[10+kNch]*adcZEMootlg[0]+corrCoeff0[10+kNch]);
645 dZEM2Corr[0] = adcZEM[1] - (corrCoeff1[11]*adcZEMoot[1]+corrCoeff0[11]);
646 dZEM2Corr[1] = adcZEMlg[1] - (corrCoeff1[11+kNch]*adcZEMootlg[1]+corrCoeff0[11+kNch]);
647 // *************************************************************************
650 if(fRecoMode==1) // p-p data
651 ReconstructEventpp(clustersTree, tZN1Corr, tZP1Corr, tZN2Corr, tZP2Corr,
652 dZEM1Corr, dZEM2Corr, sPMRef1, sPMRef2,
653 isScalerOn, scalerData, tdcData,
654 evQualityBlock, triggerBlock, chBlock, puBits);
655 else if(fRecoMode==2) // Pb-Pb data
656 ReconstructEventPbPb(clustersTree, tZN1Corr, tZP1Corr, tZN2Corr, tZP2Corr,
657 dZEM1Corr, dZEM2Corr, sPMRef1, sPMRef2,
658 isScalerOn, scalerData, tdcData,
659 evQualityBlock, triggerBlock, chBlock, puBits);
662 //_____________________________________________________________________________
663 void AliZDCReconstructor::ReconstructEventpp(TTree *clustersTree,
664 const Float_t* const corrADCZN1, const Float_t* const corrADCZP1,
665 const Float_t* const corrADCZN2, const Float_t* const corrADCZP2,
666 const Float_t* const corrADCZEM1, const Float_t* const corrADCZEM2,
667 Float_t* sPMRef1, Float_t* sPMRef2, Bool_t isScalerOn, UInt_t* scaler,
668 Int_t tdcData[32][4], const Int_t* const evQualityBlock,
669 const Int_t* const triggerBlock, const Int_t* const chBlock, UInt_t puBits) const
671 // ****************** Reconstruct one event ******************
674 /*printf("\n*************************************************\n");
675 printf(" ReconstructEventpp -> values after pedestal subtraction:\n");
676 printf(" ADCZN1 [%1.2f %1.2f %1.2f %1.2f %1.2f]\n",
677 corrADCZN1[0],corrADCZN1[1],corrADCZN1[2],corrADCZN1[3],corrADCZN1[4]);
678 printf(" ADCZP1 [%1.2f %1.2f %1.2f %1.2f %1.2f]\n",
679 corrADCZP1[0],corrADCZP1[1],corrADCZP1[2],corrADCZP1[3],corrADCZP1[4]);
680 printf(" ADCZN2 [%1.2f %1.2f %1.2f %1.2f %1.2f]\n",
681 corrADCZN2[0],corrADCZN2[1],corrADCZN2[2],corrADCZN2[3],corrADCZN2[4]);
682 printf(" ADCZP2 [%1.2f %1.2f %1.2f %1.2f %1.2f]\n",
683 corrADCZP2[0],corrADCZP2[1],corrADCZP2[2],corrADCZP2[3],corrADCZP2[4]);
684 printf(" ADCZEM1 [%1.2f] ADCZEM2 [%1.2f] \n",corrADCZEM1[0],corrADCZEM2[0]);
685 printf("*************************************************\n");*/
687 // ---------------------- Setting reco flags for ESD
689 for(Int_t ifl=0; ifl<32; ifl++) rFlags[ifl]=0;
691 if(evQualityBlock[0] == 1) rFlags[31] = 0x0;
692 else rFlags[31] = 0x1;
694 if(evQualityBlock[1] == 1) rFlags[30] = 0x1;
695 if(evQualityBlock[2] == 1) rFlags[29] = 0x1;
696 if(evQualityBlock[3] == 1) rFlags[28] = 0x1;
698 if(triggerBlock[0] == 1) rFlags[27] = 0x1;
699 if(triggerBlock[1] == 1) rFlags[26] = 0x1;
700 if(triggerBlock[2] == 1) rFlags[25] = 0x1;
701 if(triggerBlock[3] == 1) rFlags[24] = 0x1;
703 if(chBlock[0] == 1) rFlags[18] = 0x1;
704 if(chBlock[1] == 1) rFlags[17] = 0x1;
705 if(chBlock[2] == 1) rFlags[16] = 0x1;
708 rFlags[13] = puBits & 0x00000020;
709 rFlags[12] = puBits & 0x00000010;
710 rFlags[11] = puBits & 0x00000080;
711 rFlags[10] = puBits & 0x00000040;
712 rFlags[9] = puBits & 0x00000020;
713 rFlags[8] = puBits & 0x00000010;
715 if(corrADCZP1[0]>fSignalThreshold) rFlags[5] = 0x1;
716 if(corrADCZN1[0]>fSignalThreshold) rFlags[4] = 0x1;
717 if(corrADCZEM2[0]>fSignalThreshold) rFlags[3] = 0x1;
718 if(corrADCZEM1[0]>fSignalThreshold) rFlags[2] = 0x1;
719 if(corrADCZP2[0]>fSignalThreshold) rFlags[1] = 0x1;
720 if(corrADCZN2[0]>fSignalThreshold) rFlags[0] = 0x1;
722 UInt_t recoFlag = rFlags[31] << 31 | rFlags[30] << 30 | rFlags[29] << 29 | rFlags[28] << 28 |
723 rFlags[27] << 27 | rFlags[26] << 26 | rFlags[25] << 25 | rFlags[24] << 24 |
724 0x0 << 23 | 0x0 << 22 | 0x0 << 21 | 0x0 << 20 |
725 0x0 << 19 | rFlags[18] << 18 | rFlags[17] << 17 | rFlags[16] << 16 |
726 0x0 << 15 | 0x0 << 14 | rFlags[13] << 13 | rFlags[12] << 12 |
727 rFlags[11] << 11 |rFlags[10] << 10 | rFlags[9] << 9 | rFlags[8] << 8 |
728 0x0 << 7 | 0x0 << 6 | rFlags[5] << 5 | rFlags[4] << 4 |
729 rFlags[3] << 3 | rFlags[2] << 2 | rFlags[1] << 1 | rFlags[0];
730 // --------------------------------------------------
732 // ****** Retrieving calibration data
733 // --- Equalization coefficients ---------------------------------------------
734 Float_t equalCoeffZN1[5], equalCoeffZP1[5], equalCoeffZN2[5], equalCoeffZP2[5];
735 for(Int_t ji=0; ji<5; ji++){
736 equalCoeffZN1[ji] = fTowCalibData->GetZN1EqualCoeff(ji);
737 equalCoeffZP1[ji] = fTowCalibData->GetZP1EqualCoeff(ji);
738 equalCoeffZN2[ji] = fTowCalibData->GetZN2EqualCoeff(ji);
739 equalCoeffZP2[ji] = fTowCalibData->GetZP2EqualCoeff(ji);
741 // --- Energy calibration factors ------------------------------------
743 // **** Energy calibration coefficient set to 1
744 // **** (no trivial way to calibrate in p-p runs)
745 for(Int_t ij=0; ij<6; ij++) calibEne[ij] = fEnCalibData->GetEnCalib(ij);
747 // ****** Equalization of detector responses
748 Float_t equalTowZN1[10], equalTowZN2[10], equalTowZP1[10], equalTowZP2[10];
749 for(Int_t gi=0; gi<10; gi++){
751 equalTowZN1[gi] = corrADCZN1[gi]*equalCoeffZN1[gi];
752 equalTowZP1[gi] = corrADCZP1[gi]*equalCoeffZP1[gi];
753 equalTowZN2[gi] = corrADCZN2[gi]*equalCoeffZN2[gi];
754 equalTowZP2[gi] = corrADCZP2[gi]*equalCoeffZP2[gi];
757 equalTowZN1[gi] = corrADCZN1[gi]*equalCoeffZN1[gi-5];
758 equalTowZP1[gi] = corrADCZP1[gi]*equalCoeffZP1[gi-5];
759 equalTowZN2[gi] = corrADCZN2[gi]*equalCoeffZN2[gi-5];
760 equalTowZP2[gi] = corrADCZP2[gi]*equalCoeffZP2[gi-5];
764 /*printf("\n ------------- EQUALIZATION -------------\n");
765 printf(" ADCZN1 [%1.2f %1.2f %1.2f %1.2f %1.2f]\n",
766 equalTowZN1[0],equalTowZN1[1],equalTowZN1[2],equalTowZN1[3],equalTowZN1[4]);
767 printf(" ADCZP1 [%1.2f %1.2f %1.2f %1.2f %1.2f]\n",
768 equalTowZP1[0],equalTowZP1[1],equalTowZP1[2],equalTowZP1[3],equalTowZP1[4]);
769 printf(" ADCZN2 [%1.2f %1.2f %1.2f %1.2f %1.2f]\n",
770 equalTowZN2[0],equalTowZN2[1],equalTowZN2[2],equalTowZN2[3],equalTowZN2[4]);
771 printf(" ADCZP2 [%1.2f %1.2f %1.2f %1.2f %1.2f]\n",
772 equalTowZP2[0],equalTowZP2[1],equalTowZP2[2],equalTowZP2[3],equalTowZP2[4]);
773 printf(" ----------------------------------------\n");*/
775 // ****** Summed response for hadronic calorimeter (SUMMED and then CALIBRATED!)
776 Float_t calibSumZN1[]={0,0}, calibSumZN2[]={0,0}, calibSumZP1[]={0,0}, calibSumZP2[]={0,0};
777 for(Int_t gi=0; gi<5; gi++){
778 calibSumZN1[0] += equalTowZN1[gi];
779 calibSumZP1[0] += equalTowZP1[gi];
780 calibSumZN2[0] += equalTowZN2[gi];
781 calibSumZP2[0] += equalTowZP2[gi];
783 calibSumZN1[1] += equalTowZN1[gi+5];
784 calibSumZP1[1] += equalTowZP1[gi+5];
785 calibSumZN2[1] += equalTowZN2[gi+5];
786 calibSumZP2[1] += equalTowZP2[gi+5];
789 calibSumZN1[0] = calibSumZN1[0]*calibEne[0];
790 calibSumZP1[0] = calibSumZP1[0]*calibEne[1];
791 calibSumZN2[0] = calibSumZN2[0]*calibEne[2];
792 calibSumZP2[0] = calibSumZP2[0]*calibEne[3];
794 calibSumZN1[1] = calibSumZN1[1]*calibEne[0];
795 calibSumZP1[1] = calibSumZP1[1]*calibEne[1];
796 calibSumZN2[1] = calibSumZN2[1]*calibEne[2];
797 calibSumZP2[1] = calibSumZP2[1]*calibEne[3];
799 // ****** Energy calibration of detector responses
800 Float_t calibTowZN1[10], calibTowZN2[10], calibTowZP1[10], calibTowZP2[10];
801 for(Int_t gi=0; gi<5; gi++){
803 calibTowZN1[gi] = equalTowZN1[gi]*calibEne[0];
804 calibTowZP1[gi] = equalTowZP1[gi]*calibEne[1];
805 calibTowZN2[gi] = equalTowZN2[gi]*calibEne[2];
806 calibTowZP2[gi] = equalTowZP2[gi]*calibEne[3];
808 calibTowZN1[gi+5] = equalTowZN1[gi+5]*calibEne[0];
809 calibTowZP1[gi+5] = equalTowZP1[gi+5]*calibEne[1];
810 calibTowZN2[gi+5] = equalTowZN2[gi+5]*calibEne[2];
811 calibTowZP2[gi+5] = equalTowZP2[gi+5]*calibEne[3];
814 Float_t sumZEM[]={0,0}, calibZEM1[]={0,0}, calibZEM2[]={0,0};
815 calibZEM1[0] = corrADCZEM1[0]*calibEne[4];
816 calibZEM1[1] = corrADCZEM1[1]*calibEne[4];
817 calibZEM2[0] = corrADCZEM2[0]*calibEne[5];
818 calibZEM2[1] = corrADCZEM2[1]*calibEne[5];
819 for(Int_t k=0; k<2; k++) sumZEM[k] = calibZEM1[k] + calibZEM2[k];
821 /*printf("\n ------------- CALIBRATION -------------\n");
822 printf(" ADCZN1 [%1.2f %1.2f %1.2f %1.2f %1.2f]\n",
823 calibTowZN1[0],calibTowZN1[1],calibTowZN1[2],calibTowZN1[3],calibTowZN1[4]);
824 printf(" ADCZP1 [%1.2f %1.2f %1.2f %1.2f %1.2f]\n",
825 calibTowZP1[0],calibTowZP1[1],calibTowZP1[2],calibTowZP1[3],calibTowZP1[4]);
826 printf(" ADCZN2 [%1.2f %1.2f %1.2f %1.2f %1.2f]\n",
827 calibTowZN2[0],calibTowZN2[1],calibTowZN2[2],calibTowZN2[3],calibTowZN2[4]);
828 printf(" ADCZP2 [%1.2f %1.2f %1.2f %1.2f %1.2f]\n",
829 calibTowZP2[0],calibTowZP2[1],calibTowZP2[2],calibTowZP2[3],calibTowZP2[4]);
830 printf(" ADCZEM1 [%1.2f] ADCZEM2 [%1.2f] \n",calibZEM1[0],calibZEM2[0]);
831 printf(" ----------------------------------------\n");*/
833 // ****** No. of spectator and participants nucleons
834 // Variables calculated to comply with ESD structure
835 // *** N.B. -> They have a meaning only in Pb-Pb!!!!!!!!!!!!
836 Int_t nDetSpecNLeft=0, nDetSpecPLeft=0, nDetSpecNRight=0, nDetSpecPRight=0;
837 Int_t nGenSpec=0, nGenSpecLeft=0, nGenSpecRight=0;
838 Int_t nPart=0, nPartTotLeft=0, nPartTotRight=0;
839 Double_t impPar=0., impPar1=0., impPar2=0.;
841 Bool_t energyFlag = kFALSE;
842 // create the output tree
843 AliZDCReco* reco = new AliZDCReco(calibSumZN1, calibSumZP1, calibSumZN2, calibSumZP2,
844 calibTowZN1, calibTowZP1, calibTowZN2, calibTowZP2,
845 calibZEM1, calibZEM2, sPMRef1, sPMRef2,
846 nDetSpecNLeft, nDetSpecPLeft, nDetSpecNRight, nDetSpecPRight,
847 nGenSpec, nGenSpecLeft, nGenSpecRight,
848 nPart, nPartTotLeft, nPartTotRight,
849 impPar, impPar1, impPar2,
850 recoFlag, energyFlag, isScalerOn, scaler, tdcData);
852 const Int_t kBufferSize = 4000;
853 clustersTree->Branch("ZDC", "AliZDCReco", &reco, kBufferSize);
854 // write the output tree
855 clustersTree->Fill();
859 //_____________________________________________________________________________
860 void AliZDCReconstructor::ReconstructEventPbPb(TTree *clustersTree,
861 const Float_t* const corrADCZN1, const Float_t* const corrADCZP1,
862 const Float_t* const corrADCZN2, const Float_t* const corrADCZP2,
863 const Float_t* const corrADCZEM1, const Float_t* const corrADCZEM2,
864 Float_t* sPMRef1, Float_t* sPMRef2, Bool_t isScalerOn, UInt_t* scaler,
865 Int_t tdcData[32][4], const Int_t* const evQualityBlock,
866 const Int_t* const triggerBlock, const Int_t* const chBlock, UInt_t puBits) const
868 // ****************** Reconstruct one event ******************
869 // ---------------------- Setting reco flags for ESD
871 for(Int_t ifl=0; ifl<32; ifl++) rFlags[ifl]=0;
873 if(evQualityBlock[0] == 1) rFlags[31] = 0x0;
874 else rFlags[31] = 0x1;
876 if(evQualityBlock[1] == 1) rFlags[30] = 0x1;
877 if(evQualityBlock[2] == 1) rFlags[29] = 0x1;
878 if(evQualityBlock[3] == 1) rFlags[28] = 0x1;
880 if(triggerBlock[0] == 1) rFlags[27] = 0x1;
881 if(triggerBlock[1] == 1) rFlags[26] = 0x1;
882 if(triggerBlock[2] == 1) rFlags[25] = 0x1;
883 if(triggerBlock[3] == 1) rFlags[24] = 0x1;
885 if(chBlock[0] == 1) rFlags[18] = 0x1;
886 if(chBlock[1] == 1) rFlags[17] = 0x1;
887 if(chBlock[2] == 1) rFlags[16] = 0x1;
889 rFlags[13] = puBits & 0x00000020;
890 rFlags[12] = puBits & 0x00000010;
891 rFlags[11] = puBits & 0x00000080;
892 rFlags[10] = puBits & 0x00000040;
893 rFlags[9] = puBits & 0x00000020;
894 rFlags[8] = puBits & 0x00000010;
896 if(corrADCZP1[0]>fSignalThreshold) rFlags[5] = 0x1;
897 if(corrADCZN1[0]>fSignalThreshold) rFlags[4] = 0x1;
898 if(corrADCZEM2[0]>fSignalThreshold) rFlags[3] = 0x1;
899 if(corrADCZEM1[0]>fSignalThreshold) rFlags[2] = 0x1;
900 if(corrADCZP2[0]>fSignalThreshold) rFlags[1] = 0x1;
901 if(corrADCZN2[0]>fSignalThreshold) rFlags[0] = 0x1;
903 UInt_t recoFlag = rFlags[31] << 31 | rFlags[30] << 30 | rFlags[29] << 29 | rFlags[28] << 28 |
904 rFlags[27] << 27 | rFlags[26] << 26 | rFlags[25] << 25 | rFlags[24] << 24 |
905 0x0 << 23 | 0x0 << 22 | 0x0 << 21 | 0x0 << 20 |
906 0x0 << 19 | rFlags[18] << 18 | rFlags[17] << 17 | rFlags[16] << 16 |
907 0x0 << 15 | 0x0 << 14 | rFlags[13] << 13 | rFlags[12] << 12 |
908 rFlags[11] << 11 |rFlags[10] << 10 | rFlags[9] << 9 | rFlags[8] << 8 |
909 0x0 << 7 | 0x0 << 6 | rFlags[5] << 5 | rFlags[4] << 4 |
910 rFlags[3] << 3 | rFlags[2] << 2 | rFlags[1] << 1 | rFlags[0];
911 // --------------------------------------------------
915 /* printf("\n*************************************************\n");
916 printf(" ReconstructEventPbPb -> values after pedestal subtraction:\n");
917 printf(" ADCZN1 [%1.2f %1.2f %1.2f %1.2f %1.2f]\n",
918 corrADCZN1[0],corrADCZN1[1],corrADCZN1[2],corrADCZN1[3],corrADCZN1[4]);
919 printf(" ADCZP1 [%1.2f %1.2f %1.2f %1.2f %1.2f]\n",
920 corrADCZP1[0],corrADCZP1[1],corrADCZP1[2],corrADCZP1[3],corrADCZP1[4]);
921 printf(" ADCZN2 [%1.2f %1.2f %1.2f %1.2f %1.2f]\n",
922 corrADCZN2[0],corrADCZN2[1],corrADCZN2[2],corrADCZN2[3],corrADCZN2[4]);
923 printf(" ADCZP2 [%1.2f %1.2f %1.2f %1.2f %1.2f]\n",
924 corrADCZP2[0],corrADCZP2[1],corrADCZP2[2],corrADCZP2[3],corrADCZP2[4]);
925 printf(" ADCZEM1 [%1.2f] ADCZEM2 [%1.2f] \n",corrADCZEM1[0],corrADCZEM2[0]);
926 printf("*************************************************\n");
928 // ****** Retrieving calibration data
929 // --- Equalization coefficients ---------------------------------------------
930 Float_t equalCoeffZN1[5], equalCoeffZP1[5], equalCoeffZN2[5], equalCoeffZP2[5];
931 for(Int_t ji=0; ji<5; ji++){
932 equalCoeffZN1[ji] = fTowCalibData->GetZN1EqualCoeff(ji);
933 equalCoeffZP1[ji] = fTowCalibData->GetZP1EqualCoeff(ji);
934 equalCoeffZN2[ji] = fTowCalibData->GetZN2EqualCoeff(ji);
935 equalCoeffZP2[ji] = fTowCalibData->GetZP2EqualCoeff(ji);
937 // --- Energy calibration factors ------------------------------------
939 // The energy calibration object already takes into account of E_beam
940 // -> the value from the OCDB can be directly used (Jul 2010)
941 for(Int_t ij=0; ij<6; ij++) calibEne[ij] = fEnCalibData->GetEnCalib(ij);
943 // ****** Equalization of detector responses
944 Float_t equalTowZN1[10], equalTowZN2[10], equalTowZP1[10], equalTowZP2[10];
945 for(Int_t gi=0; gi<10; gi++){
947 equalTowZN1[gi] = corrADCZN1[gi]*equalCoeffZN1[gi];
948 equalTowZP1[gi] = corrADCZP1[gi]*equalCoeffZP1[gi];
949 equalTowZN2[gi] = corrADCZN2[gi]*equalCoeffZN2[gi];
950 equalTowZP2[gi] = corrADCZP2[gi]*equalCoeffZP2[gi];
953 equalTowZN1[gi] = corrADCZN1[gi]*equalCoeffZN1[gi-5];
954 equalTowZP1[gi] = corrADCZP1[gi]*equalCoeffZP1[gi-5];
955 equalTowZN2[gi] = corrADCZN2[gi]*equalCoeffZN2[gi-5];
956 equalTowZP2[gi] = corrADCZP2[gi]*equalCoeffZP2[gi-5];
961 /* printf("\n ------------- EQUALIZATION -------------\n");
962 printf(" ADCZN1 [%1.2f %1.2f %1.2f %1.2f %1.2f]\n",
963 equalTowZN1[0],equalTowZN1[1],equalTowZN1[2],equalTowZN1[3],equalTowZN1[4]);
964 printf(" ADCZP1 [%1.2f %1.2f %1.2f %1.2f %1.2f]\n",
965 equalTowZP1[0],equalTowZP1[1],equalTowZP1[2],equalTowZP1[3],equalTowZP1[4]);
966 printf(" ADCZN2 [%1.2f %1.2f %1.2f %1.2f %1.2f]\n",
967 equalTowZN2[0],equalTowZN2[1],equalTowZN2[2],equalTowZN2[3],equalTowZN2[4]);
968 printf(" ADCZP2 [%1.2f %1.2f %1.2f %1.2f %1.2f]\n",
969 equalTowZP2[0],equalTowZP2[1],equalTowZP2[2],equalTowZP2[3],equalTowZP2[4]);
970 printf(" ----------------------------------------\n");
973 // ****** Summed response for hadronic calorimeter (SUMMED and then CALIBRATED!)
974 Float_t calibSumZN1[]={0,0}, calibSumZN2[]={0,0}, calibSumZP1[]={0,0}, calibSumZP2[]={0,0};
975 for(Int_t gi=0; gi<5; gi++){
976 calibSumZN1[0] += equalTowZN1[gi];
977 calibSumZP1[0] += equalTowZP1[gi];
978 calibSumZN2[0] += equalTowZN2[gi];
979 calibSumZP2[0] += equalTowZP2[gi];
981 calibSumZN1[1] += equalTowZN1[gi+5];
982 calibSumZP1[1] += equalTowZP1[gi+5];
983 calibSumZN2[1] += equalTowZN2[gi+5];
984 calibSumZP2[1] += equalTowZP2[gi+5];
987 //fEnCalibData->Print("");
990 calibSumZN1[0] = calibSumZN1[0]*calibEne[0]*8.;
991 calibSumZP1[0] = calibSumZP1[0]*calibEne[1]*8.;
992 calibSumZN2[0] = calibSumZN2[0]*calibEne[2]*8.;
993 calibSumZP2[0] = calibSumZP2[0]*calibEne[3]*8.;
995 calibSumZN1[1] = calibSumZN1[1]*calibEne[0];
996 calibSumZP1[1] = calibSumZP1[1]*calibEne[1];
997 calibSumZN2[1] = calibSumZN2[1]*calibEne[2];
998 calibSumZP2[1] = calibSumZP2[1]*calibEne[3];
1000 Float_t sumZEM[]={0,0}, calibZEM1[]={0,0}, calibZEM2[]={0,0};
1001 calibZEM1[0] = corrADCZEM1[0]*calibEne[4]*8.;
1002 calibZEM1[1] = corrADCZEM1[1]*calibEne[4];
1003 calibZEM2[0] = corrADCZEM2[0]*calibEne[5]*8.;
1004 calibZEM2[1] = corrADCZEM2[1]*calibEne[5];
1005 for(Int_t k=0; k<2; k++) sumZEM[k] = calibZEM1[k] + calibZEM2[k];
1007 // ****** Energy calibration of detector responses
1008 Float_t calibTowZN1[10], calibTowZN2[10], calibTowZP1[10], calibTowZP2[10];
1009 for(Int_t gi=0; gi<5; gi++){
1011 calibTowZN1[gi] = equalTowZN1[gi]*2*calibEne[0]*8.;
1012 calibTowZP1[gi] = equalTowZP1[gi]*2*calibEne[1]*8.;
1013 calibTowZN2[gi] = equalTowZN2[gi]*2*calibEne[2]*8.;
1014 calibTowZP2[gi] = equalTowZP2[gi]*2*calibEne[3]*8.;
1016 calibTowZN1[gi+5] = equalTowZN1[gi+5]*2*calibEne[0];
1017 calibTowZP1[gi+5] = equalTowZP1[gi+5]*2*calibEne[1];
1018 calibTowZN2[gi+5] = equalTowZN2[gi+5]*2*calibEne[2];
1019 calibTowZP2[gi+5] = equalTowZP2[gi+5]*2*calibEne[3];
1023 /* printf("\n ------------- CALIBRATION -------------\n");
1024 printf(" ADCZN1 [%1.2f %1.2f %1.2f %1.2f %1.2f]\n",
1025 calibTowZN1[0],calibTowZN1[1],calibTowZN1[2],calibTowZN1[3],calibTowZN1[4]);
1026 printf(" ADCZP1 [%1.2f %1.2f %1.2f %1.2f %1.2f]\n",
1027 calibTowZP1[0],calibTowZP1[1],calibTowZP1[2],calibTowZP1[3],calibTowZP1[4]);
1028 printf(" ADCZN2 [%1.2f %1.2f %1.2f %1.2f %1.2f]\n",
1029 calibTowZN2[0],calibTowZN2[1],calibTowZN2[2],calibTowZN2[3],calibTowZN2[4]);
1030 printf(" ADCZP2 [%1.2f %1.2f %1.2f %1.2f %1.2f]\n",
1031 calibTowZP2[0],calibTowZP2[1],calibTowZP2[2],calibTowZP2[3],calibTowZP2[4]);
1032 printf(" ADCZEM1 [%1.2f] ADCZEM2 [%1.2f] \n",calibZEM1[0],calibZEM2[0]);
1033 printf(" ----------------------------------------\n");
1035 // ****** Number of detected spectator nucleons
1036 Int_t nDetSpecNLeft=0, nDetSpecPLeft=0, nDetSpecNRight=0, nDetSpecPRight=0;
1037 if(fBeamEnergy>0.01){
1038 nDetSpecNLeft = (Int_t) (calibSumZN1[0]/fBeamEnergy);
1039 nDetSpecPLeft = (Int_t) (calibSumZP1[0]/fBeamEnergy);
1040 nDetSpecNRight = (Int_t) (calibSumZN2[0]/fBeamEnergy);
1041 nDetSpecPRight = (Int_t) (calibSumZP2[0]/fBeamEnergy);
1043 else AliWarning(" ATTENTION!!! fBeamEnergy=0 -> N_spec will be ZERO!!! \n");
1044 /*printf("\n\t AliZDCReconstructor -> fBeamEnergy %1.0f: nDetSpecNsideA %d, nDetSpecPsideA %d,"
1045 " nDetSpecNsideC %d, nDetSpecPsideC %d\n",fBeamEnergy,nDetSpecNLeft, nDetSpecPLeft,
1046 nDetSpecNRight, nDetSpecPRight);*/
1048 Int_t nGenSpec=0, nGenSpecA=0, nGenSpecC=0;
1049 Int_t nPart=0, nPartA=0, nPartC=0;
1050 Double_t b=0., bA=0., bC=0.;
1052 if(fIsCalibrationMB == kFALSE){
1053 // ****** Reconstruction parameters ------------------
1054 if(!fgRecoParam) fgRecoParam = const_cast<AliZDCRecoParam*>(GetRecoParam());
1056 AliError(" RecoParam object not retrieved correctly: not reconstructing event!!!");
1059 TH1D* hNpartDist = fgRecoParam->GethNpartDist();
1060 TH1D* hbDist = fgRecoParam->GethbDist();
1061 Float_t fClkCenter = fgRecoParam->GetClkCenter();
1062 if(!hNpartDist || !hbDist){
1063 AliError("Something wrong in Glauber MC histos got from AliZDCREcoParamPbPb: NO EVENT RECO FOR ZDC DATA!!!\n\n");
1067 if(!fgMBCalibData) fgMBCalibData = const_cast<AliZDCMBCalib*>(GetMBCalibData());
1068 TH2F *hZDCvsZEM = fgMBCalibData->GethZDCvsZEM();
1069 TH2F *hZDCCvsZEM = fgMBCalibData->GethZDCCvsZEM();
1070 TH2F *hZDCAvsZEM = fgMBCalibData->GethZDCAvsZEM();
1072 Double_t xHighEdge = hZDCvsZEM->GetXaxis()->GetXmax();
1073 Double_t origin = xHighEdge*fClkCenter;
1075 //printf("\n\n xHighEdge %1.2f, origin %1.4f \n", xHighEdge, origin);
1077 // ====> Summed ZDC info (sideA+side C)
1078 TF1 *line = new TF1("line","[0]*x+[1]",0.,xHighEdge);
1079 Float_t y = (calibSumZN1[0]+calibSumZP1[0]+calibSumZN2[0]+calibSumZP2[0])/1000.;
1080 Float_t x = (calibZEM1[0]+calibZEM2[0])/1000.;
1081 line->SetParameter(0, y/(x-origin));
1082 line->SetParameter(1, -origin*y/(x-origin));
1084 //printf(" ***************** Summed ZDC info (sideA+side C) \n");
1085 //printf(" E_{ZEM} %1.4f, E_{ZDC} %1.2f, TF1: %1.2f*x + %1.2f ", x, y,y/(x-origin),-origin*y/(x-origin));
1087 Double_t countPerc=0;
1088 Double_t xBinCenter=0, yBinCenter=0;
1089 for(Int_t nbinx=1; nbinx<=hZDCvsZEM->GetNbinsX(); nbinx++){
1090 for(Int_t nbiny=1; nbiny<=hZDCvsZEM->GetNbinsY(); nbiny++){
1091 xBinCenter = hZDCvsZEM->GetXaxis()->GetBinCenter(nbinx);
1092 yBinCenter = hZDCvsZEM->GetYaxis()->GetBinCenter(nbiny);
1094 if(line->GetParameter(0)>0){
1095 if(yBinCenter < (line->GetParameter(0)*xBinCenter + line->GetParameter(1))){
1096 countPerc += hZDCvsZEM->GetBinContent(nbinx,nbiny);
1098 //printf(" xBinCenter %1.3f, yBinCenter %1.0f, countPerc %1.0f\n",
1099 //xBinCenter, yBinCenter, countPerc);
1103 if(yBinCenter > (line->GetParameter(0)*xBinCenter + line->GetParameter(1))){
1104 countPerc += hZDCvsZEM->GetBinContent(nbinx,nbiny);
1106 //printf(" xBinCenter %1.3f, yBinCenter %1.0f, countPerc %1.0f\n",
1107 //xBinCenter, yBinCenter, countPerc);
1113 Double_t xSecPerc = 0.;
1114 if(hZDCvsZEM->GetEntries()!=0){
1115 xSecPerc = countPerc/hZDCvsZEM->GetEntries();
1118 AliWarning(" Histogram hZDCvsZEM from OCDB has no entries!!!");
1121 //printf(" xSecPerc %1.4f \n", xSecPerc);
1124 TF1 *lineC = new TF1("lineC","[0]*x+[1]",0.,xHighEdge);
1125 Float_t yC = (calibSumZN1[0]+calibSumZP1[0])/1000.;
1126 lineC->SetParameter(0, yC/(x-origin));
1127 lineC->SetParameter(1, -origin*yC/(x-origin));
1129 //printf(" ***************** Side C \n");
1130 //printf(" E_{ZEM} %1.4f, E_{ZDCC} %1.2f, TF1: %1.2f*x + %1.2f ", x, yC,yC/(x-origin),-origin*yC/(x-origin));
1132 Double_t countPercC=0;
1133 Double_t xBinCenterC=0, yBinCenterC=0;
1134 for(Int_t nbinx=1; nbinx<=hZDCCvsZEM->GetNbinsX(); nbinx++){
1135 for(Int_t nbiny=1; nbiny<=hZDCCvsZEM->GetNbinsY(); nbiny++){
1136 xBinCenterC = hZDCCvsZEM->GetXaxis()->GetBinCenter(nbinx);
1137 yBinCenterC = hZDCCvsZEM->GetYaxis()->GetBinCenter(nbiny);
1138 if(lineC->GetParameter(0)>0){
1139 if(yBinCenterC < (lineC->GetParameter(0)*xBinCenterC + lineC->GetParameter(1))){
1140 countPercC += hZDCCvsZEM->GetBinContent(nbinx,nbiny);
1144 if(yBinCenterC > (lineC->GetParameter(0)*xBinCenterC + lineC->GetParameter(1))){
1145 countPercC += hZDCCvsZEM->GetBinContent(nbinx,nbiny);
1151 Double_t xSecPercC = 0.;
1152 if(hZDCCvsZEM->GetEntries()!=0){
1153 xSecPercC = countPercC/hZDCCvsZEM->GetEntries();
1156 AliWarning(" Histogram hZDCCvsZEM from OCDB has no entries!!!");
1159 //printf(" xSecPercC %1.4f \n", xSecPercC);
1162 TF1 *lineA = new TF1("lineA","[0]*x+[1]",0.,xHighEdge);
1163 Float_t yA = (calibSumZN2[0]+calibSumZP2[0])/1000.;
1164 lineA->SetParameter(0, yA/(x-origin));
1165 lineA->SetParameter(1, -origin*yA/(x-origin));
1168 //printf(" ***************** Side A \n");
1169 //printf(" E_{ZEM} %1.4f, E_{ZDCA} %1.2f, TF1: %1.2f*x + %1.2f ", x, yA,yA/(x-origin),-origin*yA/(x-origin));
1171 Double_t countPercA=0;
1172 Double_t xBinCenterA=0, yBinCenterA=0;
1173 for(Int_t nbinx=1; nbinx<=hZDCAvsZEM->GetNbinsX(); nbinx++){
1174 for(Int_t nbiny=1; nbiny<=hZDCAvsZEM->GetNbinsY(); nbiny++){
1175 xBinCenterA = hZDCAvsZEM->GetXaxis()->GetBinCenter(nbinx);
1176 yBinCenterA = hZDCAvsZEM->GetYaxis()->GetBinCenter(nbiny);
1177 if(lineA->GetParameter(0)>0){
1178 if(yBinCenterA < (lineA->GetParameter(0)*xBinCenterA + lineA->GetParameter(1))){
1179 countPercA += hZDCAvsZEM->GetBinContent(nbinx,nbiny);
1183 if(yBinCenterA > (lineA->GetParameter(0)*xBinCenterA + lineA->GetParameter(1))){
1184 countPercA += hZDCAvsZEM->GetBinContent(nbinx,nbiny);
1190 Double_t xSecPercA = 0.;
1191 if(hZDCAvsZEM->GetEntries()!=0){
1192 xSecPercA = countPercA/hZDCAvsZEM->GetEntries();
1195 AliWarning(" Histogram hZDCAvsZEM from OCDB has no entries!!!");
1198 //printf(" xSecPercA %1.4f \n", xSecPercA);
1200 // ****** Number of participants (from E_ZDC vs. E_ZEM correlation)
1201 Double_t nPartFrac=0., nPartFracC=0., nPartFracA=0.;
1202 for(Int_t npbin=1; npbin<hNpartDist->GetNbinsX(); npbin++){
1203 nPartFrac += (hNpartDist->GetBinContent(npbin))/(hNpartDist->GetEntries());
1204 if((1.-nPartFrac) < xSecPerc){
1205 nPart = (Int_t) hNpartDist->GetBinLowEdge(npbin);
1207 //printf(" ***************** Summed ZDC info (sideA+side C) \n");
1208 //printf(" nPartFrac %1.4f, nPart %d\n", nPartFrac, nPart);
1212 if(nPart<0) nPart=0;
1214 for(Int_t npbin=1; npbin<hNpartDist->GetNbinsX(); npbin++){
1215 nPartFracC += (hNpartDist->GetBinContent(npbin))/(hNpartDist->GetEntries());
1216 if((1.-nPartFracC) < xSecPercC){
1217 nPartC = (Int_t) hNpartDist->GetBinLowEdge(npbin);
1219 //printf(" ***************** Side C \n");
1220 //printf(" nPartFracC %1.4f, nPartC %d\n", nPartFracC, nPartC);
1224 if(nPartC<0) nPartC=0;
1226 for(Int_t npbin=1; npbin<hNpartDist->GetNbinsX(); npbin++){
1227 nPartFracA += (hNpartDist->GetBinContent(npbin))/(hNpartDist->GetEntries());
1228 if((1.-nPartFracA) < xSecPercA){
1229 nPartA = (Int_t) hNpartDist->GetBinLowEdge(npbin);
1231 //printf(" ***************** Side A \n");
1232 //printf(" nPartFracA %1.4f, nPartA %d\n\n", nPartFracA, nPartA);
1236 if(nPartA<0) nPartA=0;
1238 // ****** Impact parameter (from E_ZDC vs. E_ZEM correlation)
1239 Double_t bFrac=0., bFracC=0., bFracA=0.;
1240 for(Int_t ibbin=1; ibbin<hbDist->GetNbinsX(); ibbin++){
1241 bFrac += (hbDist->GetBinContent(ibbin))/(hbDist->GetEntries());
1242 if(bFrac > xSecPerc){
1243 b = hbDist->GetBinLowEdge(ibbin);
1248 for(Int_t ibbin=1; ibbin<hbDist->GetNbinsX(); ibbin++){
1249 bFracC += (hbDist->GetBinContent(ibbin))/(hbDist->GetEntries());
1250 if(bFracC > xSecPercC){
1251 bC = hbDist->GetBinLowEdge(ibbin);
1256 for(Int_t ibbin=1; ibbin<hbDist->GetNbinsX(); ibbin++){
1257 bFracA += (hbDist->GetBinContent(ibbin))/(hbDist->GetEntries());
1258 if(bFracA > xSecPercA){
1259 bA = hbDist->GetBinLowEdge(ibbin);
1264 // ****** Number of spectator nucleons
1265 nGenSpec = 416 - nPart;
1266 nGenSpecC = 416 - nPartC;
1267 nGenSpecA = 416 - nPartA;
1268 if(nGenSpec>416) nGenSpec=416; if(nGenSpec<0) nGenSpec=0;
1269 if(nGenSpecC>416) nGenSpecC=416; if(nGenSpecC<0) nGenSpecC=0;
1270 if(nGenSpecA>416) nGenSpecA=416; if(nGenSpecA<0) nGenSpecA=0;
1273 delete lineC; delete lineA;
1275 } // ONLY IF fIsCalibrationMB==kFALSE
1277 Bool_t energyFlag = kTRUE;
1278 AliZDCReco* reco = new AliZDCReco(calibSumZN1, calibSumZP1, calibSumZN2, calibSumZP2,
1279 calibTowZN1, calibTowZP1, calibTowZN2, calibTowZP2,
1280 calibZEM1, calibZEM2, sPMRef1, sPMRef2,
1281 nDetSpecNLeft, nDetSpecPLeft, nDetSpecNRight, nDetSpecPRight,
1282 nGenSpec, nGenSpecA, nGenSpecC,
1283 nPart, nPartA, nPartC, b, bA, bC,
1284 recoFlag, energyFlag, isScalerOn, scaler, tdcData);
1286 const Int_t kBufferSize = 4000;
1287 clustersTree->Branch("ZDC", "AliZDCReco", &reco, kBufferSize);
1289 // write the output tree
1290 clustersTree->Fill();
1295 //_____________________________________________________________________________
1296 void AliZDCReconstructor::FillZDCintoESD(TTree *clustersTree, AliESDEvent* esd) const
1298 // fill energies and number of participants to the ESD
1301 AliZDCReco* preco = &reco;
1302 clustersTree->SetBranchAddress("ZDC", &preco);
1303 clustersTree->GetEntry(0);
1305 Float_t tZN1Ene[5], tZN2Ene[5], tZP1Ene[5], tZP2Ene[5];
1306 Float_t tZN1EneLR[5], tZN2EneLR[5], tZP1EneLR[5], tZP2EneLR[5];
1307 for(Int_t i=0; i<5; i++){
1308 tZN1Ene[i] = reco.GetZN1HREnTow(i);
1309 tZN2Ene[i] = reco.GetZN2HREnTow(i);
1310 tZP1Ene[i] = reco.GetZP1HREnTow(i);
1311 tZP2Ene[i] = reco.GetZP2HREnTow(i);
1313 tZN1EneLR[i] = reco.GetZN1LREnTow(i);
1314 tZN2EneLR[i] = reco.GetZN2LREnTow(i);
1315 tZP1EneLR[i] = reco.GetZP1LREnTow(i);
1316 tZP2EneLR[i] = reco.GetZP2LREnTow(i);
1319 fESDZDC->SetZN1TowerEnergy(tZN1Ene);
1320 fESDZDC->SetZN2TowerEnergy(tZN2Ene);
1321 fESDZDC->SetZP1TowerEnergy(tZP1Ene);
1322 fESDZDC->SetZP2TowerEnergy(tZP2Ene);
1324 fESDZDC->SetZN1TowerEnergyLR(tZN1EneLR);
1325 fESDZDC->SetZN2TowerEnergyLR(tZN2EneLR);
1326 fESDZDC->SetZP1TowerEnergyLR(tZP1EneLR);
1327 fESDZDC->SetZP2TowerEnergyLR(tZP2EneLR);
1329 Int_t nPart = reco.GetNParticipants();
1330 Int_t nPartA = reco.GetNPartSideA();
1331 Int_t nPartC = reco.GetNPartSideC();
1332 Double_t b = reco.GetImpParameter();
1333 Double_t bA = reco.GetImpParSideA();
1334 Double_t bC = reco.GetImpParSideC();
1335 UInt_t recoFlag = reco.GetRecoFlag();
1337 fESDZDC->SetZDC(reco.GetZN1HREnergy(), reco.GetZP1HREnergy(),
1338 reco.GetZEM1HRsignal(), reco.GetZEM2HRsignal(),
1339 reco.GetZN2HREnergy(), reco.GetZP2HREnergy(),
1340 nPart, nPartA, nPartC, b, bA, bC, recoFlag);
1342 // Writing ZDC scaler for cross section calculation
1343 // ONLY IF the scaler has been read during the event
1344 if(reco.IsScalerOn()==kTRUE){
1346 for(Int_t jk=0; jk<32; jk++) counts[jk] = reco.GetZDCScaler(jk);
1347 fESDZDC->SetZDCScaler(counts);
1350 // Writing TDC data into ZDC ESDs
1351 Int_t tdcValues[32][4];
1352 Float_t tdcCorrected[32][4];
1353 for(Int_t jk=0; jk<32; jk++){
1354 for(Int_t lk=0; lk<4; lk++){
1355 tdcValues[jk][lk] = reco.GetZDCTDCData(jk, lk);
1358 // 4/2/2011 -> Subtracting L0 (tdcValues[15]) instead of ADC gate
1359 // we try to keep the TDC oscillations as low as possible!
1360 for(Int_t jk=0; jk<32; jk++){
1361 for(Int_t lk=0; lk<4; lk++){
1362 if(tdcValues[jk][lk]!=0.) tdcCorrected[jk][lk] = 0.025*(tdcValues[jk][lk]-tdcValues[15][0])+fMeanPhase;
1365 fESDZDC->SetZDCTDCData(tdcValues);
1366 fESDZDC->SetZDCTDCCorrected(tdcCorrected);
1367 fESDZDC->AliESDZDC::SetBit(AliESDZDC::kCorrectedTDCFilled, reco.GetEnergyFlag());
1368 fESDZDC->AliESDZDC::SetBit(AliESDZDC::kEnergyCalibratedSignal, kTRUE);
1370 if(esd) esd->SetZDCData(fESDZDC);
1373 //_____________________________________________________________________________
1374 AliCDBStorage* AliZDCReconstructor::SetStorage(const char *uri)
1376 // Setting the storage
1378 Bool_t deleteManager = kFALSE;
1380 AliCDBManager *manager = AliCDBManager::Instance();
1381 AliCDBStorage *defstorage = manager->GetDefaultStorage();
1383 if(!defstorage || !(defstorage->Contains("ZDC"))){
1384 AliWarning("No default storage set or default storage doesn't contain ZDC!");
1385 manager->SetDefaultStorage(uri);
1386 deleteManager = kTRUE;
1389 AliCDBStorage *storage = manager->GetDefaultStorage();
1392 AliCDBManager::Instance()->UnsetDefaultStorage();
1393 defstorage = 0; // the storage is killed by AliCDBManager::Instance()->Destroy()
1399 //_____________________________________________________________________________
1400 AliZDCPedestals* AliZDCReconstructor::GetPedestalData() const
1403 // Getting pedestal calibration object for ZDC set
1405 AliCDBEntry *entry = AliCDBManager::Instance()->Get("ZDC/Calib/Pedestals");
1406 if(!entry) AliFatal("No calibration data loaded!");
1407 entry->SetOwner(kFALSE);
1409 AliZDCPedestals *calibdata = dynamic_cast<AliZDCPedestals*> (entry->GetObject());
1410 if(!calibdata) AliFatal("Wrong calibration object in calibration file!");
1415 //_____________________________________________________________________________
1416 AliZDCEnCalib* AliZDCReconstructor::GetEnergyCalibData() const
1419 // Getting energy and equalization calibration object for ZDC set
1421 AliCDBEntry *entry = AliCDBManager::Instance()->Get("ZDC/Calib/EnergyCalib");
1422 if(!entry) AliFatal("No calibration data loaded!");
1423 entry->SetOwner(kFALSE);
1425 AliZDCEnCalib *calibdata = dynamic_cast<AliZDCEnCalib*> (entry->GetObject());
1426 if(!calibdata) AliFatal("Wrong calibration object in calibration file!");
1431 //_____________________________________________________________________________
1432 AliZDCTowerCalib* AliZDCReconstructor::GetTowerCalibData() const
1435 // Getting energy and equalization calibration object for ZDC set
1437 AliCDBEntry *entry = AliCDBManager::Instance()->Get("ZDC/Calib/TowerCalib");
1438 if(!entry) AliFatal("No calibration data loaded!");
1439 entry->SetOwner(kFALSE);
1441 AliZDCTowerCalib *calibdata = dynamic_cast<AliZDCTowerCalib*> (entry->GetObject());
1442 if(!calibdata) AliFatal("Wrong calibration object in calibration file!");
1447 //_____________________________________________________________________________
1448 AliZDCMBCalib* AliZDCReconstructor::GetMBCalibData() const
1451 // Getting energy and equalization calibration object for ZDC set
1453 AliCDBEntry *entry = AliCDBManager::Instance()->Get("ZDC/Calib/MBCalib");
1454 if(!entry) AliFatal("No calibration data loaded!");
1455 entry->SetOwner(kFALSE);
1457 AliZDCMBCalib *calibdata = dynamic_cast<AliZDCMBCalib*> (entry->GetObject());
1458 if(!calibdata) AliFatal("Wrong calibration object in calibration file!");