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 "AliZDCSaturationCalib.h"
45 #include "AliZDCTowerCalib.h"
46 #include "AliZDCMBCalib.h"
47 #include "AliZDCTDCCalib.h"
48 #include "AliZDCRecoParam.h"
49 #include "AliZDCRecoParampp.h"
50 #include "AliZDCRecoParamPbPb.h"
51 #include "AliRunInfo.h"
52 #include "AliLHCClockPhase.h"
55 ClassImp(AliZDCReconstructor)
56 AliZDCRecoParam *AliZDCReconstructor::fgRecoParam=0; //reconstruction parameters
57 AliZDCMBCalib *AliZDCReconstructor::fgMBCalibData=0; //calibration parameters for A-A reconstruction
59 //_____________________________________________________________________________
60 AliZDCReconstructor:: AliZDCReconstructor() :
61 fPedData(GetPedestalData()),
62 fEnCalibData(GetEnergyCalibData()),
63 fSatCalibData(GetSaturationCalibData()),
64 fTowCalibData(GetTowerCalibData()),
65 fTDCCalibData(GetTDCCalibData()),
69 fIsCalibrationMB(kFALSE),
74 // **** Default constructor
78 //_____________________________________________________________________________
79 AliZDCReconstructor::~AliZDCReconstructor()
82 // if(fgRecoParam) delete fgRecoParam;
83 if(fPedData) delete fPedData;
84 if(fEnCalibData) delete fEnCalibData;
85 if(fSatCalibData) delete fSatCalibData;
86 if(fTowCalibData) delete fTowCalibData;
87 if(fgMBCalibData) delete fgMBCalibData;
88 if(fESDZDC) delete fESDZDC;
91 //____________________________________________________________________________
92 void AliZDCReconstructor::Init()
94 // Setting reconstruction parameters
96 TString runType = GetRunInfo()->GetRunType();
97 if((runType.CompareTo("CALIBRATION_MB")) == 0){
98 fIsCalibrationMB = kTRUE;
101 TString beamType = GetRunInfo()->GetBeamType();
102 // This is a temporary solution to allow reconstruction in tests without beam
103 if(((beamType.CompareTo("UNKNOWN"))==0) &&
104 ((runType.CompareTo("PHYSICS"))==0 || (runType.CompareTo("CALIBRATION_BC"))==0)){
107 /*else if((beamType.CompareTo("UNKNOWN"))==0){
108 AliError("\t UNKNOWN beam type\n");
112 fBeamEnergy = GetRunInfo()->GetBeamEnergy();
113 if(fBeamEnergy<0.01){
114 AliWarning(" Beam energy value missing -> setting it to 1380 GeV ");
118 if(((beamType.CompareTo("pp"))==0) || ((beamType.CompareTo("p-p"))==0)
119 ||((beamType.CompareTo("PP"))==0) || ((beamType.CompareTo("P-P"))==0)){
122 else if(((beamType.CompareTo("p-A"))==0) || ((beamType.CompareTo("A-p"))==0)
123 ||((beamType.CompareTo("P-A"))==0) || ((beamType.CompareTo("A-P"))==0)){
126 else if((beamType.CompareTo("A-A")) == 0 || (beamType.CompareTo("AA")) == 0){
128 if(!fgRecoParam) fgRecoParam = const_cast<AliZDCRecoParam*>(GetRecoParam());
130 fgRecoParam->SetGlauberMCDist(fBeamEnergy);
134 AliCDBEntry *entry = AliCDBManager::Instance()->Get("GRP/Calib/LHCClockPhase");
135 if (!entry) AliFatal("LHC clock-phase shift is not found in OCDB !");
136 AliLHCClockPhase *phaseLHC = (AliLHCClockPhase*)entry->GetObject();
137 // 4/2/2011 According to A. Di Mauro BEAM1 measurement is more reliable
138 // than BEAM2 and therefore also than the average of the 2
139 fMeanPhase = phaseLHC->GetMeanPhaseB1();
141 if(fIsCalibrationMB==kFALSE)
142 AliInfo(Form("\n\n ***** ZDC reconstruction initialized for %s @ %1.0f + %1.0f GeV *****\n\n",
143 beamType.Data(), fBeamEnergy, fBeamEnergy));
145 // if EMD calibration run NO ENERGY CALIBRATION should be performed
146 // pp-like reconstruction must be performed (E cailb. coeff. = 1)
147 if((runType.CompareTo("CALIBRATION_EMD")) == 0){
152 AliInfo(Form("\n ZDC reconstruction mode %d (1 -> p-p, 2-> A-A)\n\n",fRecoMode));
154 fESDZDC = new AliESDZDC();
159 //____________________________________________________________________________
160 void AliZDCReconstructor::Init(TString beamType, Float_t beamEnergy)
162 // Setting reconstruction mode
163 // Needed to work in the HLT framework
165 fIsCalibrationMB = kFALSE;
167 fBeamEnergy = beamEnergy;
169 if(((beamType.CompareTo("pp"))==0) || ((beamType.CompareTo("p-p"))==0)
170 ||((beamType.CompareTo("PP"))==0) || ((beamType.CompareTo("P-P"))==0)){
173 else if(((beamType.CompareTo("p-A"))==0) || ((beamType.CompareTo("A-p"))==0)
174 ||((beamType.CompareTo("P-A"))==0) || ((beamType.CompareTo("A-P"))==0)){
177 else if((beamType.CompareTo("A-A")) == 0 || (beamType.CompareTo("AA")) == 0){
179 if(!fgRecoParam) fgRecoParam = const_cast<AliZDCRecoParam*>(GetRecoParam());
180 if( fgRecoParam ) fgRecoParam->SetGlauberMCDist(fBeamEnergy);
183 AliCDBEntry *entry = AliCDBManager::Instance()->Get("GRP/Calib/LHCClockPhase");
184 if (!entry) AliFatal("LHC clock-phase shift is not found in OCDB !");
185 AliLHCClockPhase *phaseLHC = (AliLHCClockPhase*)entry->GetObject();
186 fMeanPhase = phaseLHC->GetMeanPhase();
188 fESDZDC = new AliESDZDC();
190 AliInfo(Form("\n\n ***** ZDC reconstruction initialized for %s @ %1.0f + %1.0f GeV *****\n\n",
191 beamType.Data(), fBeamEnergy, fBeamEnergy));
195 //_____________________________________________________________________________
196 void AliZDCReconstructor::Reconstruct(TTree* digitsTree, TTree* clustersTree) const
198 // *** Local ZDC reconstruction for digits
199 // Works on the current event
201 // Retrieving calibration data
202 // Parameters for mean value pedestal subtraction
204 Float_t meanPed[2*kNch];
205 for(Int_t jj=0; jj<2*kNch; jj++) meanPed[jj] = fPedData->GetMeanPed(jj);
206 // Parameters pedestal subtraction through correlation with out-of-time signals
207 Float_t corrCoeff0[2*kNch], corrCoeff1[2*kNch];
208 for(Int_t jj=0; jj<2*kNch; jj++){
209 corrCoeff0[jj] = fPedData->GetPedCorrCoeff0(jj);
210 corrCoeff1[jj] = fPedData->GetPedCorrCoeff1(jj);
215 AliZDCDigit* pdigit = &digit;
216 digitsTree->SetBranchAddress("ZDC", &pdigit);
217 //printf("\n\t # of digits in tree: %d\n",(Int_t) digitsTree->GetEntries());
220 Float_t tZN1Corr[10], tZP1Corr[10], tZN2Corr[10], tZP2Corr[10];
221 Float_t dZEM1Corr[2], dZEM2Corr[2], sPMRef1[2], sPMRef2[2];
222 for(Int_t i=0; i<10; i++){
223 tZN1Corr[i] = tZP1Corr[i] = tZN2Corr[i] = tZP2Corr[i] = 0.;
224 if(i<2) dZEM1Corr[i] = dZEM2Corr[i] = sPMRef1[i] = sPMRef2[i] = 0.;
227 Int_t digNentries = digitsTree->GetEntries();
228 Float_t ootDigi[kNch]; Int_t i=0;
229 // -- Reading out-of-time signals (last kNch entries) for current event
231 for(Int_t iDigit=kNch; iDigit<digNentries; iDigit++){
232 if(i<=kNch) ootDigi[i] = digitsTree->GetEntry(iDigit);
233 else AliWarning(" Can't read more out of time values: index>kNch !!!\n");
238 for(Int_t iDigit=0; iDigit<(digNentries/2); iDigit++) {
239 digitsTree->GetEntry(iDigit);
240 if (!pdigit) continue;
242 Int_t det = digit.GetSector(0);
243 Int_t quad = digit.GetSector(1);
245 Float_t ped2SubHg=0., ped2SubLg=0.;
247 if(det==1) pedindex = quad;
248 else if(det==2) pedindex = quad+5;
249 else if(det==3) pedindex = quad+9;
250 else if(det==4) pedindex = quad+12;
251 else if(det==5) pedindex = quad+17;
253 else pedindex = (det-1)/3+22;
256 ped2SubHg = meanPed[pedindex];
257 ped2SubLg = meanPed[pedindex+kNch];
259 else if(fPedSubMode==1){
260 ped2SubHg = corrCoeff1[pedindex]*ootDigi[pedindex]+corrCoeff0[pedindex];
261 ped2SubLg = corrCoeff1[pedindex+kNch]*ootDigi[pedindex+kNch]+corrCoeff0[pedindex+kNch];
264 if(quad != 5){ // ZDC (not reference PTMs!)
265 if(det == 1){ // *** ZNC
266 tZN1Corr[quad] = (Float_t) (digit.GetADCValue(0)-ped2SubHg);
267 tZN1Corr[quad+5] = (Float_t) (digit.GetADCValue(1)-ped2SubLg);
269 else if(det == 2){ // *** ZP1
270 tZP1Corr[quad] = (Float_t) (digit.GetADCValue(0)-ped2SubHg);
271 tZP1Corr[quad+5] = (Float_t) (digit.GetADCValue(1)-ped2SubLg);
274 if(quad == 1){ // *** ZEM1
275 dZEM1Corr[0] += (Float_t) (digit.GetADCValue(0)-ped2SubHg);
276 dZEM1Corr[1] += (Float_t) (digit.GetADCValue(1)-ped2SubLg);
278 else if(quad == 2){ // *** ZEM2
279 dZEM2Corr[0] += (Float_t) (digit.GetADCValue(0)-ped2SubHg);
280 dZEM2Corr[1] += (Float_t) (digit.GetADCValue(1)-ped2SubLg);
283 else if(det == 4){ // *** ZN2
284 tZN2Corr[quad] = (Float_t) (digit.GetADCValue(0)-ped2SubHg);
285 tZN2Corr[quad+5] = (Float_t) (digit.GetADCValue(1)-ped2SubLg);
287 else if(det == 5){ // *** ZP2
288 tZP2Corr[quad] = (Float_t) (digit.GetADCValue(0)-ped2SubHg);
289 tZP2Corr[quad+5] = (Float_t) (digit.GetADCValue(1)-ped2SubLg);
292 else{ // Reference PMs
294 sPMRef1[0] = (Float_t) (digit.GetADCValue(0)-ped2SubHg);
295 sPMRef1[1] = (Float_t) (digit.GetADCValue(1)-ped2SubLg);
298 sPMRef2[0] = (Float_t) (digit.GetADCValue(0)-ped2SubHg);
299 sPMRef2[1] = (Float_t) (digit.GetADCValue(1)-ped2SubLg);
304 /*printf("AliZDCReconstructor: digit #%d det %d quad %d pedHG %1.0f pedLG %1.0f\n",
305 iDigit, det, quad, ped2SubHg, ped2SubLg);
306 printf(" -> pedindex %d\n", pedindex);
307 printf(" HGChain -> RawDig %d DigCorr %1.2f",
308 digit.GetADCValue(0), digit.GetADCValue(0)-ped2SubHg);
309 printf(" LGChain -> RawDig %d DigCorr %1.2f\n",
310 digit.GetADCValue(1), digit.GetADCValue(1)-ped2SubLg);*/
316 for(Int_t jj=0; jj<32; jj++){
318 for(Int_t ii=0; ii<4; ii++) tdc[jj][ii]=0;
321 Int_t evQualityBlock[4] = {1,0,0,0};
322 Int_t triggerBlock[4] = {0,0,0,0};
323 Int_t chBlock[3] = {0,0,0};
326 // reconstruct the event
328 ReconstructEventpp(clustersTree, tZN1Corr, tZP1Corr, tZN2Corr, tZP2Corr,
329 dZEM1Corr, dZEM2Corr, sPMRef1, sPMRef2,
331 evQualityBlock, triggerBlock, chBlock, puBits);
332 else if(fRecoMode==2)
333 ReconstructEventPbPb(clustersTree, tZN1Corr, tZP1Corr, tZN2Corr, tZP2Corr,
334 dZEM1Corr, dZEM2Corr, sPMRef1, sPMRef2,
336 evQualityBlock, triggerBlock, chBlock, puBits);
339 //_____________________________________________________________________________
340 void AliZDCReconstructor::Reconstruct(AliRawReader* rawReader, TTree* clustersTree) const
342 // *** ZDC raw data reconstruction
343 // Works on the current event
345 // Retrieving calibration data
346 // Parameters for pedestal subtraction
348 Float_t meanPed[2*kNch];
349 for(Int_t jj=0; jj<2*kNch; jj++) meanPed[jj] = fPedData->GetMeanPed(jj);
350 // Parameters pedestal subtraction through correlation with out-of-time signals
351 Float_t corrCoeff0[2*kNch], corrCoeff1[2*kNch];
352 for(Int_t jj=0; jj<2*kNch; jj++){
353 corrCoeff0[jj] = fPedData->GetPedCorrCoeff0(jj);
354 corrCoeff1[jj] = fPedData->GetPedCorrCoeff1(jj);
355 //printf(" %d %1.4f %1.4f\n", jj,corrCoeff0[jj],corrCoeff1[jj]);
358 Int_t adcZN1[5], adcZN1oot[5], adcZN1lg[5], adcZN1ootlg[5];
359 Int_t adcZP1[5], adcZP1oot[5], adcZP1lg[5], adcZP1ootlg[5];
360 Int_t adcZN2[5], adcZN2oot[5], adcZN2lg[5], adcZN2ootlg[5];
361 Int_t adcZP2[5], adcZP2oot[5], adcZP2lg[5], adcZP2ootlg[5];
362 Int_t adcZEM[2], adcZEMoot[2], adcZEMlg[2], adcZEMootlg[2];
363 Int_t pmRef[2], pmRefoot[2], pmReflg[2], pmRefootlg[2];
364 for(Int_t ich=0; ich<5; ich++){
365 adcZN1[ich] = adcZN1oot[ich] = adcZN1lg[ich] = adcZN1ootlg[ich] = 0;
366 adcZP1[ich] = adcZP1oot[ich] = adcZP1lg[ich] = adcZP1ootlg[ich] = 0;
367 adcZN2[ich] = adcZN2oot[ich] = adcZN2lg[ich] = adcZN2ootlg[ich] = 0;
368 adcZP2[ich] = adcZP2oot[ich] = adcZP2lg[ich] = adcZP2ootlg[ich] = 0;
370 adcZEM[ich] = adcZEMoot[ich] = adcZEMlg[ich] = adcZEMootlg[ich] = 0;
371 pmRef[ich] = pmRefoot[ich] = pmReflg[ich] = pmRefootlg[ich] = 0;
375 Float_t tZN1Corr[10], tZP1Corr[10], tZN2Corr[10], tZP2Corr[10];
376 Float_t dZEM1Corr[2], dZEM2Corr[2], sPMRef1[2], sPMRef2[2];
377 for(Int_t i=0; i<10; i++){
378 tZN1Corr[i] = tZP1Corr[i] = tZN2Corr[i] = tZP2Corr[i] = 0.;
379 if(i<2) dZEM1Corr[i] = dZEM2Corr[i] = sPMRef1[i] = sPMRef2[i] = 0.;
382 Bool_t isScalerOn=kFALSE;
383 Int_t jsc=0, itdc=0, iprevtdc=-1, ihittdc=0;
384 UInt_t scalerData[32];
385 Int_t tdcData[32][4];
386 for(Int_t k=0; k<32; k++){
388 for(Int_t i=0; i<4; i++) tdcData[k][i]=0;
392 Int_t evQualityBlock[4] = {1,0,0,0};
393 Int_t triggerBlock[4] = {0,0,0,0};
394 Int_t chBlock[3] = {0,0,0};
397 Int_t kFirstADCGeo=0, kLastADCGeo=3, kScalerGeo=8, kZDCTDCGeo=4, kPUGeo=29;
398 //Int_t kTrigScales=30, kTrigHistory=31;
400 // loop over raw data
401 //rawReader->Reset();
402 AliZDCRawStream rawData(rawReader);
403 while(rawData.Next()){
405 // ***************************** Reading ADCs
406 if((rawData.GetADCModule()>=kFirstADCGeo) && (rawData.GetADCModule()<=kLastADCGeo)){
407 //printf(" **** Reading ADC raw data from module %d **** \n",rawData.GetADCModule());
409 if((rawData.IsADCDataWord()) && (rawData.GetNChannelsOn()<48)) chBlock[0] = kTRUE;
410 if((rawData.IsADCDataWord()) && (rawData.IsOverflow() == kTRUE)) chBlock[1] = kTRUE;
411 if((rawData.IsADCDataWord()) && (rawData.IsUnderflow() == kTRUE)) chBlock[2] = kTRUE;
412 if((rawData.IsADCDataWord()) && (rawData.IsADCEventGood() == kTRUE)) evQualityBlock[0] = kTRUE;
414 if((rawData.IsADCDataWord()) && (rawData.IsUnderflow()==kFALSE)
415 && (rawData.IsOverflow()==kFALSE) && (rawData.IsADCEventGood()==kTRUE)){
417 Int_t adcMod = rawData.GetADCModule();
418 Int_t det = rawData.GetSector(0);
419 Int_t quad = rawData.GetSector(1);
420 Int_t gain = rawData.GetADCGain();
423 // Mean pedestal value subtraction -------------------------------------------------------
424 if(fPedSubMode == 0){
425 // **** Pb-Pb data taking 2010 -> subtracting some ch. from correlation ****
426 // Not interested in o.o.t. signals (ADC modules 2, 3)
427 //if(adcMod == 2 || adcMod == 3) continue;
428 // **** Pb-Pb data taking 2011 -> subtracting only ZEM from correlation ****
430 if(adcMod==0 || adcMod==1){
431 if(gain==0) adcZEM[quad-1] = rawData.GetADCValue();
432 else adcZEMlg[quad-1] = rawData.GetADCValue();
434 else if(adcMod==2 || adcMod==3){
435 if(gain==0) adcZEMoot[quad-1] = rawData.GetADCValue();
436 else adcZEMootlg[quad-1] = rawData.GetADCValue();
439 // When oot values are read the ADC modules 2, 3 can be skipped!!!
440 if(adcMod == 2 || adcMod == 3) continue;
442 // *************************************************************************
443 if(quad != 5){ // ZDCs (not reference PTMs)
446 if(gain == 0) tZN1Corr[quad] += (Float_t) (rawData.GetADCValue()-meanPed[pedindex]);
447 else tZN1Corr[quad+5] += (Float_t) (rawData.GetADCValue()-meanPed[pedindex+kNch]);
451 if(gain == 0) tZP1Corr[quad] += (Float_t) (rawData.GetADCValue()-meanPed[pedindex]);
452 else tZP1Corr[quad+5] += (Float_t) (rawData.GetADCValue()-meanPed[pedindex+kNch]);
457 if(gain == 0) dZEM1Corr[0] += (Float_t) (rawData.GetADCValue()-meanPed[pedindex]);
458 else dZEM1Corr[1] += (Float_t) (rawData.GetADCValue()-meanPed[pedindex+kNch]);
461 if(gain == 0) dZEM2Corr[0] += (Float_t) (rawData.GetADCValue()-meanPed[pedindex]);
462 else dZEM2Corr[1] += (Float_t) (rawData.GetADCValue()-meanPed[pedindex+kNch]);
467 if(gain == 0) tZN2Corr[quad] += (Float_t) (rawData.GetADCValue()-meanPed[pedindex]);
468 else tZN2Corr[quad+5] += (Float_t) (rawData.GetADCValue()-meanPed[pedindex+kNch]);
472 if(gain == 0) tZP2Corr[quad] += (Float_t) (rawData.GetADCValue()-meanPed[pedindex]);
473 else tZP2Corr[quad+5] += (Float_t) (rawData.GetADCValue()-meanPed[pedindex+kNch]);
476 else{ // reference PM
477 pedindex = (det-1)/3 + 22;
479 if(gain==0) sPMRef1[0] += (Float_t) (rawData.GetADCValue()-meanPed[pedindex]);
480 else sPMRef1[1] += (Float_t) (rawData.GetADCValue()-meanPed[pedindex+kNch]);
483 if(gain==0) sPMRef2[0] += (Float_t) (rawData.GetADCValue()-meanPed[pedindex]);
484 else sPMRef2[1] += (Float_t) (rawData.GetADCValue()-meanPed[pedindex+kNch]);
489 printf(" AliZDCReconstructor: det %d quad %d res %d -> Pedestal[%d] %1.0f",
490 det,quad,gain, pedindex, meanPed[pedindex]);
491 printf(" RawADC %d ADCCorr %1.0f\n",
492 rawData.GetADCValue(), rawData.GetADCValue()-meanPed[pedindex]);
494 }// mean pedestal subtraction
495 // Pedestal subtraction from correlation ------------------------------------------------
496 else if(fPedSubMode == 1){
498 if(adcMod==0 || adcMod==1){
499 if(quad != 5){ // signals from ZDCs
501 if(gain==0) adcZN1[quad] = rawData.GetADCValue();
502 else adcZN1lg[quad] = rawData.GetADCValue();
505 if(gain==0) adcZP1[quad] = rawData.GetADCValue();
506 else adcZP1lg[quad] = rawData.GetADCValue();
509 if(gain==0) adcZEM[quad-1] = rawData.GetADCValue();
510 else adcZEMlg[quad-1] = rawData.GetADCValue();
513 if(gain==0) adcZN2[quad] = rawData.GetADCValue();
514 else adcZN2lg[quad] = rawData.GetADCValue();
517 if(gain==0) adcZP2[quad] = rawData.GetADCValue();
518 else adcZP2lg[quad] = rawData.GetADCValue();
521 else{ // signals from reference PM
522 if(gain==0) pmRef[quad-1] = rawData.GetADCValue();
523 else pmReflg[quad-1] = rawData.GetADCValue();
526 // Out-of-time pedestals
527 else if(adcMod==2 || adcMod==3){
528 if(quad != 5){ // signals from ZDCs
530 if(gain==0) adcZN1oot[quad] = rawData.GetADCValue();
531 else adcZN1ootlg[quad] = rawData.GetADCValue();
534 if(gain==0) adcZP1oot[quad] = rawData.GetADCValue();
535 else adcZP1ootlg[quad] = rawData.GetADCValue();
538 if(gain==0) adcZEMoot[quad-1] = rawData.GetADCValue();
539 else adcZEMootlg[quad-1] = rawData.GetADCValue();
542 if(gain==0) adcZN2oot[quad] = rawData.GetADCValue();
543 else adcZN2ootlg[quad] = rawData.GetADCValue();
546 if(gain==0) adcZP2oot[quad] = rawData.GetADCValue();
547 else adcZP2ootlg[quad] = rawData.GetADCValue();
550 else{ // signals from reference PM
551 if(gain==0) pmRefoot[quad-1] = rawData.GetADCValue();
552 else pmRefootlg[quad-1] = rawData.GetADCValue();
555 } // pedestal subtraction from correlation
557 /*printf("\t AliZDCReconstructor: det %d quad %d res %d -> Ped[%d] = %1.0f\n",
558 det,quad,gain, pedindex, meanPed[pedindex]);*/
561 // ***************************** Reading Scaler
562 else if(rawData.GetADCModule()==kScalerGeo){
563 if(rawData.IsScalerWord()==kTRUE){
565 scalerData[jsc] = rawData.GetTriggerCount();
567 //printf(" Reconstructed VME Scaler: %d %d ",jsc,scalerData[jsc]);
572 // ***************************** Reading ZDC TDC
573 else if(rawData.GetADCModule()==kZDCTDCGeo && rawData.IsZDCTDCDatum()==kTRUE){
574 itdc = rawData.GetChannel();
575 if(itdc==iprevtdc) ihittdc++;
578 if(ihittdc<4) tdcData[itdc][ihittdc] = rawData.GetZDCTDCDatum();
580 //if(ihittdc==0) printf(" TDC%d %d ",itdc, tdcData[itdc][ihittdc]);
582 // ***************************** Reading PU
583 else if(rawData.GetADCModule()==kPUGeo){
584 puBits = rawData.GetDetectorPattern();
586 // ***************************** Reading trigger history
587 else if(rawData.IstriggerHistoryWord()==kTRUE){
588 triggerBlock[0] = rawData.IsCPTInputEMDTrigger();
589 triggerBlock[1] = rawData.IsCPTInputSemiCentralTrigger();
590 triggerBlock[2] = rawData.IsCPTInputCentralTrigger();
591 triggerBlock[3] = rawData.IsCPTInputMBTrigger();
597 for(Int_t t=0; t<5; t++){
598 tZN1Corr[t] = adcZN1[t] - (corrCoeff1[t]*adcZN1oot[t]+corrCoeff0[t]);
599 tZN1Corr[t+5] = adcZN1lg[t] - (corrCoeff1[t+kNch]*adcZN1ootlg[t]+corrCoeff0[t+kNch]);
601 tZP1Corr[t] = adcZP1[t] - (corrCoeff1[t+5]*adcZP1oot[t]+corrCoeff0[t+5]);
602 tZP1Corr[t+5] = adcZP1lg[t] - (corrCoeff1[t+5+kNch]*adcZP1ootlg[t]+corrCoeff0[t+5+kNch]);
604 tZN2Corr[t] = adcZN2[t] - (corrCoeff1[t+12]*adcZN2oot[t]+corrCoeff0[t+12]);
605 tZN2Corr[t+5] = adcZN2lg[t] - (corrCoeff1[t+12+kNch]*adcZN2ootlg[t]+corrCoeff0[t+12+kNch]);
607 tZP2Corr[t] = adcZP2[t] - (corrCoeff1[t+17]*adcZP2oot[t]+corrCoeff0[t+17]);
608 tZP2Corr[t+5] = adcZP2lg[t] - (corrCoeff1[t+17+kNch]*adcZP2ootlg[t]+corrCoeff0[t+17+kNch]);
610 dZEM1Corr[0] = adcZEM[0] - (corrCoeff1[10]*adcZEMoot[0]+corrCoeff0[10]);
611 dZEM1Corr[1] = adcZEMlg[0] - (corrCoeff1[10+kNch]*adcZEMootlg[0]+corrCoeff0[10+kNch]);
612 dZEM2Corr[0] = adcZEM[1] - (corrCoeff1[11]*adcZEMoot[1]+corrCoeff0[11]);
613 dZEM2Corr[1] = adcZEMlg[1] - (corrCoeff1[11+kNch]*adcZEMootlg[1]+corrCoeff0[11+kNch]);
615 sPMRef1[0] = pmRef[0] - (corrCoeff1[22]*pmRefoot[0]+corrCoeff0[22]);
616 sPMRef1[1] = pmReflg[0] - (corrCoeff1[22+kNch]*pmRefootlg[0]+corrCoeff0[22+kNch]);
617 sPMRef2[0] = pmRef[0] - (corrCoeff1[23]*pmRefoot[1]+corrCoeff0[23]);
618 sPMRef2[1] = pmReflg[0] - (corrCoeff1[23+kNch]*pmRefootlg[1]+corrCoeff0[23+kNch]);
620 if(fPedSubMode==0 && fRecoMode==2){
621 // **** Pb-Pb data taking 2011 -> subtracting some ch. from correlation ****
622 //tZN1Corr[0] = adcZN1[0] - (corrCoeff1[0]*adcZN1oot[0]+corrCoeff0[0]);
623 //tZN1Corr[5] = adcZN1lg[0] - (corrCoeff1[kNch]*adcZN1ootlg[0]+corrCoeff0[kNch]);
625 //printf(" adcZN1 %d adcZN1oot %d tZN1Corr %1.2f \n", adcZN1[0],adcZN1oot[0],tZN1Corr[0]);
626 //printf(" adcZN1lg %d adcZN1ootlg %d tZN1Corrlg %1.2f \n", adcZN1lg[0],adcZN1ootlg[0],tZN1Corr[5]);
628 //tZP1Corr[2] = adcZP1[2] - (corrCoeff1[2+5]*adcZP1oot[2]+corrCoeff0[2+5]);
629 //tZP1Corr[2+5] = adcZP1lg[2] - (corrCoeff1[2+5+kNch]*adcZP1ootlg[2]+corrCoeff0[2+5+kNch]);
631 dZEM1Corr[0] = adcZEM[0] - (corrCoeff1[10]*adcZEMoot[0]+corrCoeff0[10]);
632 dZEM1Corr[1] = adcZEMlg[0] - (corrCoeff1[10+kNch]*adcZEMootlg[0]+corrCoeff0[10+kNch]);
633 dZEM2Corr[0] = adcZEM[1] - (corrCoeff1[11]*adcZEMoot[1]+corrCoeff0[11]);
634 dZEM2Corr[1] = adcZEMlg[1] - (corrCoeff1[11+kNch]*adcZEMootlg[1]+corrCoeff0[11+kNch]);
635 // *************************************************************************
637 else if(fPedSubMode==0 && fRecoMode==1){
638 // **** p-p data taking 2011 -> temporary patch to overcome DA problem ****
640 dZEM1Corr[0] = adcZEM[0] - meanPed[10];
641 dZEM1Corr[1] = adcZEMlg[0] - meanPed[10+kNch];
642 dZEM2Corr[0] = adcZEM[1] - meanPed[11];
643 dZEM2Corr[1] = adcZEMlg[1] - meanPed[11+kNch];
644 // *************************************************************************
647 if(fRecoMode==1) // p-p data
648 ReconstructEventpp(clustersTree, tZN1Corr, tZP1Corr, tZN2Corr, tZP2Corr,
649 dZEM1Corr, dZEM2Corr, sPMRef1, sPMRef2,
650 isScalerOn, scalerData, tdcData,
651 evQualityBlock, triggerBlock, chBlock, puBits);
652 else if(fRecoMode==2) // Pb-Pb data
653 ReconstructEventPbPb(clustersTree, tZN1Corr, tZP1Corr, tZN2Corr, tZP2Corr,
654 dZEM1Corr, dZEM2Corr, sPMRef1, sPMRef2,
655 isScalerOn, scalerData, tdcData,
656 evQualityBlock, triggerBlock, chBlock, puBits);
659 //_____________________________________________________________________________
660 void AliZDCReconstructor::ReconstructEventpp(TTree *clustersTree,
661 const Float_t* const corrADCZN1, const Float_t* const corrADCZP1,
662 const Float_t* const corrADCZN2, const Float_t* const corrADCZP2,
663 const Float_t* const corrADCZEM1, const Float_t* const corrADCZEM2,
664 Float_t* sPMRef1, Float_t* sPMRef2, Bool_t isScalerOn, UInt_t* scaler,
665 Int_t tdcData[32][4], const Int_t* const evQualityBlock,
666 const Int_t* const triggerBlock, const Int_t* const chBlock, UInt_t puBits) const
668 // ****************** Reconstruct one event ******************
671 /*printf("\n*************************************************\n");
672 printf(" ReconstructEventpp -> values after pedestal subtraction:\n");
673 printf(" ADCZN1 [%1.2f %1.2f %1.2f %1.2f %1.2f]\n",
674 corrADCZN1[0],corrADCZN1[1],corrADCZN1[2],corrADCZN1[3],corrADCZN1[4]);
675 printf(" ADCZP1 [%1.2f %1.2f %1.2f %1.2f %1.2f]\n",
676 corrADCZP1[0],corrADCZP1[1],corrADCZP1[2],corrADCZP1[3],corrADCZP1[4]);
677 printf(" ADCZN2 [%1.2f %1.2f %1.2f %1.2f %1.2f]\n",
678 corrADCZN2[0],corrADCZN2[1],corrADCZN2[2],corrADCZN2[3],corrADCZN2[4]);
679 printf(" ADCZP2 [%1.2f %1.2f %1.2f %1.2f %1.2f]\n",
680 corrADCZP2[0],corrADCZP2[1],corrADCZP2[2],corrADCZP2[3],corrADCZP2[4]);
681 printf(" ADCZEM1 [%1.2f] ADCZEM2 [%1.2f] \n",corrADCZEM1[0],corrADCZEM2[0]);
682 printf("*************************************************\n");*/
684 // ---------------------- Setting reco flags for ESD
686 for(Int_t ifl=0; ifl<32; ifl++) rFlags[ifl]=0;
688 if(evQualityBlock[0] == 1) rFlags[31] = 0x0;
689 else rFlags[31] = 0x1;
691 if(evQualityBlock[1] == 1) rFlags[30] = 0x1;
692 if(evQualityBlock[2] == 1) rFlags[29] = 0x1;
693 if(evQualityBlock[3] == 1) rFlags[28] = 0x1;
695 if(triggerBlock[0] == 1) rFlags[27] = 0x1;
696 if(triggerBlock[1] == 1) rFlags[26] = 0x1;
697 if(triggerBlock[2] == 1) rFlags[25] = 0x1;
698 if(triggerBlock[3] == 1) rFlags[24] = 0x1;
700 if(chBlock[0] == 1) rFlags[18] = 0x1;
701 if(chBlock[1] == 1) rFlags[17] = 0x1;
702 if(chBlock[2] == 1) rFlags[16] = 0x1;
705 rFlags[13] = puBits & 0x00000020;
706 rFlags[12] = puBits & 0x00000010;
707 rFlags[11] = puBits & 0x00000080;
708 rFlags[10] = puBits & 0x00000040;
709 rFlags[9] = puBits & 0x00000020;
710 rFlags[8] = puBits & 0x00000010;
712 if(corrADCZP1[0]>fSignalThreshold) rFlags[5] = 0x1;
713 if(corrADCZN1[0]>fSignalThreshold) rFlags[4] = 0x1;
714 if(corrADCZEM2[0]>fSignalThreshold) rFlags[3] = 0x1;
715 if(corrADCZEM1[0]>fSignalThreshold) rFlags[2] = 0x1;
716 if(corrADCZP2[0]>fSignalThreshold) rFlags[1] = 0x1;
717 if(corrADCZN2[0]>fSignalThreshold) rFlags[0] = 0x1;
719 UInt_t recoFlag = rFlags[31] << 31 | rFlags[30] << 30 | rFlags[29] << 29 | rFlags[28] << 28 |
720 rFlags[27] << 27 | rFlags[26] << 26 | rFlags[25] << 25 | rFlags[24] << 24 |
721 0x0 << 23 | 0x0 << 22 | 0x0 << 21 | 0x0 << 20 |
722 0x0 << 19 | rFlags[18] << 18 | rFlags[17] << 17 | rFlags[16] << 16 |
723 0x0 << 15 | 0x0 << 14 | rFlags[13] << 13 | rFlags[12] << 12 |
724 rFlags[11] << 11 |rFlags[10] << 10 | rFlags[9] << 9 | rFlags[8] << 8 |
725 0x0 << 7 | 0x0 << 6 | rFlags[5] << 5 | rFlags[4] << 4 |
726 rFlags[3] << 3 | rFlags[2] << 2 | rFlags[1] << 1 | rFlags[0];
727 // --------------------------------------------------
729 // ****** Retrieving calibration data
730 // --- Equalization coefficients ---------------------------------------------
731 Float_t equalCoeffZN1[5], equalCoeffZP1[5], equalCoeffZN2[5], equalCoeffZP2[5];
732 for(Int_t ji=0; ji<5; ji++){
733 equalCoeffZN1[ji] = fTowCalibData->GetZN1EqualCoeff(ji);
734 equalCoeffZP1[ji] = fTowCalibData->GetZP1EqualCoeff(ji);
735 equalCoeffZN2[ji] = fTowCalibData->GetZN2EqualCoeff(ji);
736 equalCoeffZP2[ji] = fTowCalibData->GetZP2EqualCoeff(ji);
738 // --- Energy calibration factors ------------------------------------
739 Float_t calibEne[6], calibSatZNA[4], calibSatZNC[4];
740 // **** Energy calibration coefficient set to 1
741 // **** (no trivial way to calibrate in p-p runs)
742 for(Int_t ij=0; ij<6; ij++) calibEne[ij] = fEnCalibData->GetEnCalib(ij);
743 for(Int_t ij=0; ij<4; ij++){
744 calibSatZNA[ij] = fSatCalibData->GetZNASatCalib(ij);
745 calibSatZNC[ij] = fSatCalibData->GetZNCSatCalib(ij);
748 // ****** Equalization of detector responses
749 Float_t equalTowZN1[10], equalTowZN2[10], equalTowZP1[10], equalTowZP2[10];
750 for(Int_t gi=0; gi<10; gi++){
752 equalTowZN1[gi] = corrADCZN1[gi]*equalCoeffZN1[gi];
753 equalTowZP1[gi] = corrADCZP1[gi]*equalCoeffZP1[gi];
754 equalTowZN2[gi] = corrADCZN2[gi]*equalCoeffZN2[gi];
755 equalTowZP2[gi] = corrADCZP2[gi]*equalCoeffZP2[gi];
758 equalTowZN1[gi] = corrADCZN1[gi]*equalCoeffZN1[gi-5];
759 equalTowZP1[gi] = corrADCZP1[gi]*equalCoeffZP1[gi-5];
760 equalTowZN2[gi] = corrADCZN2[gi]*equalCoeffZN2[gi-5];
761 equalTowZP2[gi] = corrADCZP2[gi]*equalCoeffZP2[gi-5];
765 /*printf("\n ------------- EQUALIZATION -------------\n");
766 printf(" ADCZN1 [%1.2f %1.2f %1.2f %1.2f %1.2f]\n",
767 equalTowZN1[0],equalTowZN1[1],equalTowZN1[2],equalTowZN1[3],equalTowZN1[4]);
768 printf(" ADCZP1 [%1.2f %1.2f %1.2f %1.2f %1.2f]\n",
769 equalTowZP1[0],equalTowZP1[1],equalTowZP1[2],equalTowZP1[3],equalTowZP1[4]);
770 printf(" ADCZN2 [%1.2f %1.2f %1.2f %1.2f %1.2f]\n",
771 equalTowZN2[0],equalTowZN2[1],equalTowZN2[2],equalTowZN2[3],equalTowZN2[4]);
772 printf(" ADCZP2 [%1.2f %1.2f %1.2f %1.2f %1.2f]\n",
773 equalTowZP2[0],equalTowZP2[1],equalTowZP2[2],equalTowZP2[3],equalTowZP2[4]);
774 printf(" ----------------------------------------\n");*/
776 // *** p-A RUN 2013 -> new calibration object
777 // to take into account saturation in ZN PMC
778 // -> 5th order pol. fun. to be applied BEFORE en. calibration
779 equalTowZN1[0] = equalTowZN1[0] + calibSatZNC[0]*equalTowZN1[0]*equalTowZN1[0] +
780 calibSatZNC[1]*equalTowZN1[0]*equalTowZN1[0]*equalTowZN1[0] +
781 calibSatZNC[2]*equalTowZN1[0]*equalTowZN1[0]*equalTowZN1[0]*equalTowZN1[0] +
782 calibSatZNC[3]*equalTowZN1[0]*equalTowZN1[0]*equalTowZN1[0]*equalTowZN1[0]*equalTowZN1[0];
783 equalTowZN2[0] = equalTowZN2[0] + calibSatZNA[0]*equalTowZN2[0]*equalTowZN2[0] +
784 calibSatZNA[1]*equalTowZN2[0]*equalTowZN2[0]*equalTowZN2[0] +
785 calibSatZNA[2]*equalTowZN2[0]*equalTowZN2[0]*equalTowZN2[0]*equalTowZN2[0] +
786 calibSatZNA[3]*equalTowZN2[0]*equalTowZN2[0]*equalTowZN2[0]*equalTowZN2[0]*equalTowZN2[0];
789 /*printf("\n ------------- SATURATION CORRECTION -------------\n");
790 printf(" ZNC PMC %1.2f\n", equalTowZN1[0]);
791 printf(" ZNA PMC %1.2f\n", equalTowZN2[0]);
792 printf(" ----------------------------------------\n");*/
794 // ****** Summed response for hadronic calorimeter (SUMMED and then CALIBRATED!)
795 Float_t calibSumZN1[]={0,0}, calibSumZN2[]={0,0}, calibSumZP1[]={0,0}, calibSumZP2[]={0,0};
796 for(Int_t gi=0; gi<5; gi++){
797 calibSumZN1[0] += equalTowZN1[gi];
798 calibSumZP1[0] += equalTowZP1[gi];
799 calibSumZN2[0] += equalTowZN2[gi];
800 calibSumZP2[0] += equalTowZP2[gi];
802 calibSumZN1[1] += equalTowZN1[gi+5];
803 calibSumZP1[1] += equalTowZP1[gi+5];
804 calibSumZN2[1] += equalTowZN2[gi+5];
805 calibSumZP2[1] += equalTowZP2[gi+5];
808 calibSumZN1[0] = calibSumZN1[0]*calibEne[0];
809 calibSumZP1[0] = calibSumZP1[0]*calibEne[1];
810 calibSumZN2[0] = calibSumZN2[0]*calibEne[2];
811 calibSumZP2[0] = calibSumZP2[0]*calibEne[3];
813 calibSumZN1[1] = calibSumZN1[1]*calibEne[0];
814 calibSumZP1[1] = calibSumZP1[1]*calibEne[1];
815 calibSumZN2[1] = calibSumZN2[1]*calibEne[2];
816 calibSumZP2[1] = calibSumZP2[1]*calibEne[3];
818 // ****** Energy calibration of detector responses
819 Float_t calibTowZN1[10], calibTowZN2[10], calibTowZP1[10], calibTowZP2[10];
820 for(Int_t gi=0; gi<5; gi++){
822 calibTowZN1[gi] = equalTowZN1[gi]*calibEne[0];
823 calibTowZP1[gi] = equalTowZP1[gi]*calibEne[1];
824 calibTowZN2[gi] = equalTowZN2[gi]*calibEne[2];
825 calibTowZP2[gi] = equalTowZP2[gi]*calibEne[3];
827 calibTowZN1[gi+5] = equalTowZN1[gi+5]*calibEne[0];
828 calibTowZP1[gi+5] = equalTowZP1[gi+5]*calibEne[1];
829 calibTowZN2[gi+5] = equalTowZN2[gi+5]*calibEne[2];
830 calibTowZP2[gi+5] = equalTowZP2[gi+5]*calibEne[3];
833 Float_t calibZEM1[]={0,0}, calibZEM2[]={0,0};
834 calibZEM1[0] = corrADCZEM1[0]*calibEne[4];
835 calibZEM1[1] = corrADCZEM1[1]*calibEne[4];
836 calibZEM2[0] = corrADCZEM2[0]*calibEne[5];
837 calibZEM2[1] = corrADCZEM2[1]*calibEne[5];
839 /*printf("\n ------------- CALIBRATION -------------\n");
840 printf(" ADCZN1 [%1.2f %1.2f %1.2f %1.2f %1.2f]\n",
841 calibTowZN1[0],calibTowZN1[1],calibTowZN1[2],calibTowZN1[3],calibTowZN1[4]);
842 printf(" ADCZP1 [%1.2f %1.2f %1.2f %1.2f %1.2f]\n",
843 calibTowZP1[0],calibTowZP1[1],calibTowZP1[2],calibTowZP1[3],calibTowZP1[4]);
844 printf(" ADCZN2 [%1.2f %1.2f %1.2f %1.2f %1.2f]\n",
845 calibTowZN2[0],calibTowZN2[1],calibTowZN2[2],calibTowZN2[3],calibTowZN2[4]);
846 printf(" ADCZP2 [%1.2f %1.2f %1.2f %1.2f %1.2f]\n",
847 calibTowZP2[0],calibTowZP2[1],calibTowZP2[2],calibTowZP2[3],calibTowZP2[4]);
848 printf(" ADCZEM1 [%1.2f] ADCZEM2 [%1.2f] \n",calibZEM1[0],calibZEM2[0]);
849 printf(" ----------------------------------------\n");*/
851 // ****** No. of spectator and participants nucleons
852 // Variables calculated to comply with ESD structure
853 // *** N.B. -> They have a meaning only in Pb-Pb!!!!!!!!!!!!
854 Int_t nDetSpecNLeft=0, nDetSpecPLeft=0, nDetSpecNRight=0, nDetSpecPRight=0;
855 Int_t nGenSpec=0, nGenSpecLeft=0, nGenSpecRight=0;
856 Int_t nPart=0, nPartTotLeft=0, nPartTotRight=0;
857 Double_t impPar=0., impPar1=0., impPar2=0.;
859 Bool_t energyFlag = kFALSE;
860 // create the output tree
861 AliZDCReco* reco = new AliZDCReco(calibSumZN1, calibSumZP1, calibSumZN2, calibSumZP2,
862 calibTowZN1, calibTowZP1, calibTowZN2, calibTowZP2,
863 calibZEM1, calibZEM2, sPMRef1, sPMRef2,
864 nDetSpecNLeft, nDetSpecPLeft, nDetSpecNRight, nDetSpecPRight,
865 nGenSpec, nGenSpecLeft, nGenSpecRight,
866 nPart, nPartTotLeft, nPartTotRight,
867 impPar, impPar1, impPar2,
868 recoFlag, energyFlag, isScalerOn, scaler, tdcData);
870 const Int_t kBufferSize = 4000;
871 clustersTree->Branch("ZDC", "AliZDCReco", &reco, kBufferSize);
872 // write the output tree
873 clustersTree->Fill();
877 //_____________________________________________________________________________
878 void AliZDCReconstructor::ReconstructEventPbPb(TTree *clustersTree,
879 const Float_t* const corrADCZN1, const Float_t* const corrADCZP1,
880 const Float_t* const corrADCZN2, const Float_t* const corrADCZP2,
881 const Float_t* const corrADCZEM1, const Float_t* const corrADCZEM2,
882 Float_t* sPMRef1, Float_t* sPMRef2, Bool_t isScalerOn, UInt_t* scaler,
883 Int_t tdcData[32][4], const Int_t* const evQualityBlock,
884 const Int_t* const triggerBlock, const Int_t* const chBlock, UInt_t puBits) const
886 // ****************** Reconstruct one event ******************
887 // ---------------------- Setting reco flags for ESD
889 for(Int_t ifl=0; ifl<32; ifl++) rFlags[ifl]=0;
891 if(evQualityBlock[0] == 1) rFlags[31] = 0x0;
892 else rFlags[31] = 0x1;
894 if(evQualityBlock[1] == 1) rFlags[30] = 0x1;
895 if(evQualityBlock[2] == 1) rFlags[29] = 0x1;
896 if(evQualityBlock[3] == 1) rFlags[28] = 0x1;
898 if(triggerBlock[0] == 1) rFlags[27] = 0x1;
899 if(triggerBlock[1] == 1) rFlags[26] = 0x1;
900 if(triggerBlock[2] == 1) rFlags[25] = 0x1;
901 if(triggerBlock[3] == 1) rFlags[24] = 0x1;
903 if(chBlock[0] == 1) rFlags[18] = 0x1;
904 if(chBlock[1] == 1) rFlags[17] = 0x1;
905 if(chBlock[2] == 1) rFlags[16] = 0x1;
907 rFlags[13] = puBits & 0x00000020;
908 rFlags[12] = puBits & 0x00000010;
909 rFlags[11] = puBits & 0x00000080;
910 rFlags[10] = puBits & 0x00000040;
911 rFlags[9] = puBits & 0x00000020;
912 rFlags[8] = puBits & 0x00000010;
914 if(corrADCZP1[0]>fSignalThreshold) rFlags[5] = 0x1;
915 if(corrADCZN1[0]>fSignalThreshold) rFlags[4] = 0x1;
916 if(corrADCZEM2[0]>fSignalThreshold) rFlags[3] = 0x1;
917 if(corrADCZEM1[0]>fSignalThreshold) rFlags[2] = 0x1;
918 if(corrADCZP2[0]>fSignalThreshold) rFlags[1] = 0x1;
919 if(corrADCZN2[0]>fSignalThreshold) rFlags[0] = 0x1;
921 UInt_t recoFlag = rFlags[31] << 31 | rFlags[30] << 30 | rFlags[29] << 29 | rFlags[28] << 28 |
922 rFlags[27] << 27 | rFlags[26] << 26 | rFlags[25] << 25 | rFlags[24] << 24 |
923 0x0 << 23 | 0x0 << 22 | 0x0 << 21 | 0x0 << 20 |
924 0x0 << 19 | rFlags[18] << 18 | rFlags[17] << 17 | rFlags[16] << 16 |
925 0x0 << 15 | 0x0 << 14 | rFlags[13] << 13 | rFlags[12] << 12 |
926 rFlags[11] << 11 |rFlags[10] << 10 | rFlags[9] << 9 | rFlags[8] << 8 |
927 0x0 << 7 | 0x0 << 6 | rFlags[5] << 5 | rFlags[4] << 4 |
928 rFlags[3] << 3 | rFlags[2] << 2 | rFlags[1] << 1 | rFlags[0];
929 // --------------------------------------------------
933 /* printf("\n*************************************************\n");
934 printf(" ReconstructEventPbPb -> values after pedestal subtraction:\n");
935 printf(" ADCZN1 [%1.2f %1.2f %1.2f %1.2f %1.2f]\n",
936 corrADCZN1[0],corrADCZN1[1],corrADCZN1[2],corrADCZN1[3],corrADCZN1[4]);
937 printf(" ADCZP1 [%1.2f %1.2f %1.2f %1.2f %1.2f]\n",
938 corrADCZP1[0],corrADCZP1[1],corrADCZP1[2],corrADCZP1[3],corrADCZP1[4]);
939 printf(" ADCZN2 [%1.2f %1.2f %1.2f %1.2f %1.2f]\n",
940 corrADCZN2[0],corrADCZN2[1],corrADCZN2[2],corrADCZN2[3],corrADCZN2[4]);
941 printf(" ADCZP2 [%1.2f %1.2f %1.2f %1.2f %1.2f]\n",
942 corrADCZP2[0],corrADCZP2[1],corrADCZP2[2],corrADCZP2[3],corrADCZP2[4]);
943 printf(" ADCZEM1 [%1.2f] ADCZEM2 [%1.2f] \n",corrADCZEM1[0],corrADCZEM2[0]);
944 printf("*************************************************\n");
946 // ****** Retrieving calibration data
947 // --- Equalization coefficients ---------------------------------------------
948 Float_t equalCoeffZN1[5], equalCoeffZP1[5], equalCoeffZN2[5], equalCoeffZP2[5];
949 for(Int_t ji=0; ji<5; ji++){
950 equalCoeffZN1[ji] = fTowCalibData->GetZN1EqualCoeff(ji);
951 equalCoeffZP1[ji] = fTowCalibData->GetZP1EqualCoeff(ji);
952 equalCoeffZN2[ji] = fTowCalibData->GetZN2EqualCoeff(ji);
953 equalCoeffZP2[ji] = fTowCalibData->GetZP2EqualCoeff(ji);
955 // --- Energy calibration factors ------------------------------------
956 Float_t calibEne[6], calibSatZNA[4], calibSatZNC[4];
957 // **** Energy calibration coefficient set to 1
958 // **** (no trivial way to calibrate in p-p runs)
959 for(Int_t ij=0; ij<6; ij++) calibEne[ij] = fEnCalibData->GetEnCalib(ij);
960 for(Int_t ij=0; ij<4; ij++){
961 calibSatZNA[ij] = fSatCalibData->GetZNASatCalib(ij);
962 calibSatZNC[ij] = fSatCalibData->GetZNCSatCalib(ij);
965 // ****** Equalization of detector responses
966 Float_t equalTowZN1[10], equalTowZN2[10], equalTowZP1[10], equalTowZP2[10];
967 for(Int_t gi=0; gi<10; gi++){
969 equalTowZN1[gi] = corrADCZN1[gi]*equalCoeffZN1[gi];
970 equalTowZP1[gi] = corrADCZP1[gi]*equalCoeffZP1[gi];
971 equalTowZN2[gi] = corrADCZN2[gi]*equalCoeffZN2[gi];
972 equalTowZP2[gi] = corrADCZP2[gi]*equalCoeffZP2[gi];
975 equalTowZN1[gi] = corrADCZN1[gi]*equalCoeffZN1[gi-5];
976 equalTowZP1[gi] = corrADCZP1[gi]*equalCoeffZP1[gi-5];
977 equalTowZN2[gi] = corrADCZN2[gi]*equalCoeffZN2[gi-5];
978 equalTowZP2[gi] = corrADCZP2[gi]*equalCoeffZP2[gi-5];
983 /* printf("\n ------------- EQUALIZATION -------------\n");
984 printf(" ADCZN1 [%1.2f %1.2f %1.2f %1.2f %1.2f]\n",
985 equalTowZN1[0],equalTowZN1[1],equalTowZN1[2],equalTowZN1[3],equalTowZN1[4]);
986 printf(" ADCZP1 [%1.2f %1.2f %1.2f %1.2f %1.2f]\n",
987 equalTowZP1[0],equalTowZP1[1],equalTowZP1[2],equalTowZP1[3],equalTowZP1[4]);
988 printf(" ADCZN2 [%1.2f %1.2f %1.2f %1.2f %1.2f]\n",
989 equalTowZN2[0],equalTowZN2[1],equalTowZN2[2],equalTowZN2[3],equalTowZN2[4]);
990 printf(" ADCZP2 [%1.2f %1.2f %1.2f %1.2f %1.2f]\n",
991 equalTowZP2[0],equalTowZP2[1],equalTowZP2[2],equalTowZP2[3],equalTowZP2[4]);
992 printf(" ----------------------------------------\n");
995 // *** p-A RUN 2013 -> new calibration object
996 // to take into account saturation in ZN PMC
997 // -> 5th order pol. fun. to be applied BEFORE en. calibration
998 equalTowZN1[0] = equalTowZN1[0] + calibSatZNC[0]*equalTowZN1[0]*equalTowZN1[0] +
999 calibSatZNC[1]*equalTowZN1[0]*equalTowZN1[0]*equalTowZN1[0] +
1000 calibSatZNC[2]*equalTowZN1[0]*equalTowZN1[0]*equalTowZN1[0]*equalTowZN1[0] +
1001 calibSatZNC[3]*equalTowZN1[0]*equalTowZN1[0]*equalTowZN1[0]*equalTowZN1[0]*equalTowZN1[0];
1002 equalTowZN2[0] = equalTowZN2[0] + calibSatZNA[0]*equalTowZN2[0]*equalTowZN2[0] +
1003 calibSatZNA[1]*equalTowZN2[0]*equalTowZN2[0]*equalTowZN2[0] +
1004 calibSatZNA[2]*equalTowZN2[0]*equalTowZN2[0]*equalTowZN2[0]*equalTowZN2[0] +
1005 calibSatZNA[3]*equalTowZN2[0]*equalTowZN2[0]*equalTowZN2[0]*equalTowZN2[0]*equalTowZN2[0];
1007 // ****** Summed response for hadronic calorimeter (SUMMED and then CALIBRATED!)
1008 Float_t calibSumZN1[]={0,0}, calibSumZN2[]={0,0}, calibSumZP1[]={0,0}, calibSumZP2[]={0,0};
1009 for(Int_t gi=0; gi<5; gi++){
1010 calibSumZN1[0] += equalTowZN1[gi];
1011 calibSumZP1[0] += equalTowZP1[gi];
1012 calibSumZN2[0] += equalTowZN2[gi];
1013 calibSumZP2[0] += equalTowZP2[gi];
1015 calibSumZN1[1] += equalTowZN1[gi+5];
1016 calibSumZP1[1] += equalTowZP1[gi+5];
1017 calibSumZN2[1] += equalTowZN2[gi+5];
1018 calibSumZP2[1] += equalTowZP2[gi+5];
1021 //fEnCalibData->Print("");
1024 calibSumZN1[0] = calibSumZN1[0]*calibEne[0]*8.;
1025 calibSumZP1[0] = calibSumZP1[0]*calibEne[1]*8.;
1026 calibSumZN2[0] = calibSumZN2[0]*calibEne[2]*8.;
1027 calibSumZP2[0] = calibSumZP2[0]*calibEne[3]*8.;
1029 calibSumZN1[1] = calibSumZN1[1]*calibEne[0];
1030 calibSumZP1[1] = calibSumZP1[1]*calibEne[1];
1031 calibSumZN2[1] = calibSumZN2[1]*calibEne[2];
1032 calibSumZP2[1] = calibSumZP2[1]*calibEne[3];
1034 Float_t calibZEM1[]={0,0}, calibZEM2[]={0,0};
1035 calibZEM1[0] = corrADCZEM1[0]*calibEne[4];
1036 calibZEM1[1] = corrADCZEM1[1]*calibEne[4];
1037 calibZEM2[0] = corrADCZEM2[0]*calibEne[5];
1038 calibZEM2[1] = corrADCZEM2[1]*calibEne[5];
1040 // ****** Energy calibration of detector responses
1041 Float_t calibTowZN1[10], calibTowZN2[10], calibTowZP1[10], calibTowZP2[10];
1042 for(Int_t gi=0; gi<5; gi++){
1044 calibTowZN1[gi] = equalTowZN1[gi]*2*calibEne[0]*8.;
1045 calibTowZP1[gi] = equalTowZP1[gi]*2*calibEne[1]*8.;
1046 calibTowZN2[gi] = equalTowZN2[gi]*2*calibEne[2]*8.;
1047 calibTowZP2[gi] = equalTowZP2[gi]*2*calibEne[3]*8.;
1049 calibTowZN1[gi+5] = equalTowZN1[gi+5]*2*calibEne[0];
1050 calibTowZP1[gi+5] = equalTowZP1[gi+5]*2*calibEne[1];
1051 calibTowZN2[gi+5] = equalTowZN2[gi+5]*2*calibEne[2];
1052 calibTowZP2[gi+5] = equalTowZP2[gi+5]*2*calibEne[3];
1056 /* printf("\n ------------- CALIBRATION -------------\n");
1057 printf(" ADCZN1 [%1.2f %1.2f %1.2f %1.2f %1.2f]\n",
1058 calibTowZN1[0],calibTowZN1[1],calibTowZN1[2],calibTowZN1[3],calibTowZN1[4]);
1059 printf(" ADCZP1 [%1.2f %1.2f %1.2f %1.2f %1.2f]\n",
1060 calibTowZP1[0],calibTowZP1[1],calibTowZP1[2],calibTowZP1[3],calibTowZP1[4]);
1061 printf(" ADCZN2 [%1.2f %1.2f %1.2f %1.2f %1.2f]\n",
1062 calibTowZN2[0],calibTowZN2[1],calibTowZN2[2],calibTowZN2[3],calibTowZN2[4]);
1063 printf(" ADCZP2 [%1.2f %1.2f %1.2f %1.2f %1.2f]\n",
1064 calibTowZP2[0],calibTowZP2[1],calibTowZP2[2],calibTowZP2[3],calibTowZP2[4]);
1065 printf(" ADCZEM1 [%1.2f] ADCZEM2 [%1.2f] \n",calibZEM1[0],calibZEM2[0]);
1066 printf(" ----------------------------------------\n");
1068 // ****** Number of detected spectator nucleons
1069 Int_t nDetSpecNLeft=0, nDetSpecPLeft=0, nDetSpecNRight=0, nDetSpecPRight=0;
1070 if(fBeamEnergy>0.01){
1071 nDetSpecNLeft = (Int_t) (calibSumZN1[0]/fBeamEnergy);
1072 nDetSpecPLeft = (Int_t) (calibSumZP1[0]/fBeamEnergy);
1073 nDetSpecNRight = (Int_t) (calibSumZN2[0]/fBeamEnergy);
1074 nDetSpecPRight = (Int_t) (calibSumZP2[0]/fBeamEnergy);
1076 else AliWarning(" ATTENTION!!! fBeamEnergy=0 -> N_spec will be ZERO!!! \n");
1077 /*printf("\n\t AliZDCReconstructor -> fBeamEnergy %1.0f: nDetSpecNsideA %d, nDetSpecPsideA %d,"
1078 " nDetSpecNsideC %d, nDetSpecPsideC %d\n",fBeamEnergy,nDetSpecNLeft, nDetSpecPLeft,
1079 nDetSpecNRight, nDetSpecPRight);*/
1081 Int_t nGenSpec=0, nGenSpecA=0, nGenSpecC=0;
1082 Int_t nPart=0, nPartA=0, nPartC=0;
1083 Double_t b=0., bA=0., bC=0.;
1085 if(fIsCalibrationMB == kFALSE){
1086 // ****** Reconstruction parameters ------------------
1087 if(!fgRecoParam) fgRecoParam = const_cast<AliZDCRecoParam*>(GetRecoParam());
1089 AliError(" RecoParam object not retrieved correctly: not reconstructing ZDC event!!!");
1092 TH1D* hNpartDist = fgRecoParam->GethNpartDist();
1093 TH1D* hbDist = fgRecoParam->GethbDist();
1094 Float_t fClkCenter = fgRecoParam->GetClkCenter();
1095 if(!hNpartDist || !hbDist){
1096 AliError("Something wrong in Glauber MC histos got from AliZDCREcoParamPbPb: NO EVENT RECO FOR ZDC DATA!!!\n\n");
1100 if(!fgMBCalibData) fgMBCalibData = const_cast<AliZDCMBCalib*>(GetMBCalibData());
1101 TH2F *hZDCvsZEM = fgMBCalibData->GethZDCvsZEM();
1102 TH2F *hZDCCvsZEM = fgMBCalibData->GethZDCCvsZEM();
1103 TH2F *hZDCAvsZEM = fgMBCalibData->GethZDCAvsZEM();
1105 Double_t xHighEdge = hZDCvsZEM->GetXaxis()->GetXmax();
1106 Double_t origin = xHighEdge*fClkCenter;
1108 //printf("\n\n xHighEdge %1.2f, origin %1.4f \n", xHighEdge, origin);
1110 // ====> Summed ZDC info (sideA+side C)
1111 TF1 *line = new TF1("line","[0]*x+[1]",0.,xHighEdge);
1112 Float_t y = (calibSumZN1[0]+calibSumZP1[0]+calibSumZN2[0]+calibSumZP2[0])/1000.;
1113 Float_t x = (calibZEM1[0]+calibZEM2[0])/1000.;
1114 line->SetParameter(0, y/(x-origin));
1115 line->SetParameter(1, -origin*y/(x-origin));
1117 //printf(" ***************** Summed ZDC info (sideA+side C) \n");
1118 //printf(" E_{ZEM} %1.4f, E_{ZDC} %1.2f, TF1: %1.2f*x + %1.2f ", x, y,y/(x-origin),-origin*y/(x-origin));
1120 Double_t countPerc=0;
1121 Double_t xBinCenter=0, yBinCenter=0;
1122 for(Int_t nbinx=1; nbinx<=hZDCvsZEM->GetNbinsX(); nbinx++){
1123 for(Int_t nbiny=1; nbiny<=hZDCvsZEM->GetNbinsY(); nbiny++){
1124 xBinCenter = hZDCvsZEM->GetXaxis()->GetBinCenter(nbinx);
1125 yBinCenter = hZDCvsZEM->GetYaxis()->GetBinCenter(nbiny);
1127 if(line->GetParameter(0)>0){
1128 if(yBinCenter < (line->GetParameter(0)*xBinCenter + line->GetParameter(1))){
1129 countPerc += hZDCvsZEM->GetBinContent(nbinx,nbiny);
1131 //printf(" xBinCenter %1.3f, yBinCenter %1.0f, countPerc %1.0f\n",
1132 //xBinCenter, yBinCenter, countPerc);
1136 if(yBinCenter > (line->GetParameter(0)*xBinCenter + line->GetParameter(1))){
1137 countPerc += hZDCvsZEM->GetBinContent(nbinx,nbiny);
1139 //printf(" xBinCenter %1.3f, yBinCenter %1.0f, countPerc %1.0f\n",
1140 //xBinCenter, yBinCenter, countPerc);
1146 Double_t xSecPerc = 0.;
1147 if(hZDCvsZEM->GetEntries()!=0){
1148 xSecPerc = countPerc/hZDCvsZEM->GetEntries();
1151 AliWarning(" Histogram hZDCvsZEM from OCDB has no entries!!!");
1154 //printf(" xSecPerc %1.4f \n", xSecPerc);
1157 TF1 *lineC = new TF1("lineC","[0]*x+[1]",0.,xHighEdge);
1158 Float_t yC = (calibSumZN1[0]+calibSumZP1[0])/1000.;
1159 lineC->SetParameter(0, yC/(x-origin));
1160 lineC->SetParameter(1, -origin*yC/(x-origin));
1162 //printf(" ***************** Side C \n");
1163 //printf(" E_{ZEM} %1.4f, E_{ZDCC} %1.2f, TF1: %1.2f*x + %1.2f ", x, yC,yC/(x-origin),-origin*yC/(x-origin));
1165 Double_t countPercC=0;
1166 Double_t xBinCenterC=0, yBinCenterC=0;
1167 for(Int_t nbinx=1; nbinx<=hZDCCvsZEM->GetNbinsX(); nbinx++){
1168 for(Int_t nbiny=1; nbiny<=hZDCCvsZEM->GetNbinsY(); nbiny++){
1169 xBinCenterC = hZDCCvsZEM->GetXaxis()->GetBinCenter(nbinx);
1170 yBinCenterC = hZDCCvsZEM->GetYaxis()->GetBinCenter(nbiny);
1171 if(lineC->GetParameter(0)>0){
1172 if(yBinCenterC < (lineC->GetParameter(0)*xBinCenterC + lineC->GetParameter(1))){
1173 countPercC += hZDCCvsZEM->GetBinContent(nbinx,nbiny);
1177 if(yBinCenterC > (lineC->GetParameter(0)*xBinCenterC + lineC->GetParameter(1))){
1178 countPercC += hZDCCvsZEM->GetBinContent(nbinx,nbiny);
1184 Double_t xSecPercC = 0.;
1185 if(hZDCCvsZEM->GetEntries()!=0){
1186 xSecPercC = countPercC/hZDCCvsZEM->GetEntries();
1189 AliWarning(" Histogram hZDCCvsZEM from OCDB has no entries!!!");
1192 //printf(" xSecPercC %1.4f \n", xSecPercC);
1195 TF1 *lineA = new TF1("lineA","[0]*x+[1]",0.,xHighEdge);
1196 Float_t yA = (calibSumZN2[0]+calibSumZP2[0])/1000.;
1197 lineA->SetParameter(0, yA/(x-origin));
1198 lineA->SetParameter(1, -origin*yA/(x-origin));
1201 //printf(" ***************** Side A \n");
1202 //printf(" E_{ZEM} %1.4f, E_{ZDCA} %1.2f, TF1: %1.2f*x + %1.2f ", x, yA,yA/(x-origin),-origin*yA/(x-origin));
1204 Double_t countPercA=0;
1205 Double_t xBinCenterA=0, yBinCenterA=0;
1206 for(Int_t nbinx=1; nbinx<=hZDCAvsZEM->GetNbinsX(); nbinx++){
1207 for(Int_t nbiny=1; nbiny<=hZDCAvsZEM->GetNbinsY(); nbiny++){
1208 xBinCenterA = hZDCAvsZEM->GetXaxis()->GetBinCenter(nbinx);
1209 yBinCenterA = hZDCAvsZEM->GetYaxis()->GetBinCenter(nbiny);
1210 if(lineA->GetParameter(0)>0){
1211 if(yBinCenterA < (lineA->GetParameter(0)*xBinCenterA + lineA->GetParameter(1))){
1212 countPercA += hZDCAvsZEM->GetBinContent(nbinx,nbiny);
1216 if(yBinCenterA > (lineA->GetParameter(0)*xBinCenterA + lineA->GetParameter(1))){
1217 countPercA += hZDCAvsZEM->GetBinContent(nbinx,nbiny);
1223 Double_t xSecPercA = 0.;
1224 if(hZDCAvsZEM->GetEntries()!=0){
1225 xSecPercA = countPercA/hZDCAvsZEM->GetEntries();
1228 AliWarning(" Histogram hZDCAvsZEM from OCDB has no entries!!!");
1231 //printf(" xSecPercA %1.4f \n", xSecPercA);
1233 // ****** Number of participants (from E_ZDC vs. E_ZEM correlation)
1234 Double_t nPartFrac=0., nPartFracC=0., nPartFracA=0.;
1235 for(Int_t npbin=1; npbin<hNpartDist->GetNbinsX(); npbin++){
1236 nPartFrac += (hNpartDist->GetBinContent(npbin))/(hNpartDist->GetEntries());
1237 if((1.-nPartFrac) < xSecPerc){
1238 nPart = (Int_t) hNpartDist->GetBinLowEdge(npbin);
1240 //printf(" ***************** Summed ZDC info (sideA+side C) \n");
1241 //printf(" nPartFrac %1.4f, nPart %d\n", nPartFrac, nPart);
1245 if(nPart<0) nPart=0;
1247 for(Int_t npbin=1; npbin<hNpartDist->GetNbinsX(); npbin++){
1248 nPartFracC += (hNpartDist->GetBinContent(npbin))/(hNpartDist->GetEntries());
1249 if((1.-nPartFracC) < xSecPercC){
1250 nPartC = (Int_t) hNpartDist->GetBinLowEdge(npbin);
1252 //printf(" ***************** Side C \n");
1253 //printf(" nPartFracC %1.4f, nPartC %d\n", nPartFracC, nPartC);
1257 if(nPartC<0) nPartC=0;
1259 for(Int_t npbin=1; npbin<hNpartDist->GetNbinsX(); npbin++){
1260 nPartFracA += (hNpartDist->GetBinContent(npbin))/(hNpartDist->GetEntries());
1261 if((1.-nPartFracA) < xSecPercA){
1262 nPartA = (Int_t) hNpartDist->GetBinLowEdge(npbin);
1264 //printf(" ***************** Side A \n");
1265 //printf(" nPartFracA %1.4f, nPartA %d\n\n", nPartFracA, nPartA);
1269 if(nPartA<0) nPartA=0;
1271 // ****** Impact parameter (from E_ZDC vs. E_ZEM correlation)
1272 Double_t bFrac=0., bFracC=0., bFracA=0.;
1273 for(Int_t ibbin=1; ibbin<hbDist->GetNbinsX(); ibbin++){
1274 bFrac += (hbDist->GetBinContent(ibbin))/(hbDist->GetEntries());
1275 if(bFrac > xSecPerc){
1276 b = hbDist->GetBinLowEdge(ibbin);
1281 for(Int_t ibbin=1; ibbin<hbDist->GetNbinsX(); ibbin++){
1282 bFracC += (hbDist->GetBinContent(ibbin))/(hbDist->GetEntries());
1283 if(bFracC > xSecPercC){
1284 bC = hbDist->GetBinLowEdge(ibbin);
1289 for(Int_t ibbin=1; ibbin<hbDist->GetNbinsX(); ibbin++){
1290 bFracA += (hbDist->GetBinContent(ibbin))/(hbDist->GetEntries());
1291 if(bFracA > xSecPercA){
1292 bA = hbDist->GetBinLowEdge(ibbin);
1297 // ****** Number of spectator nucleons
1298 nGenSpec = 416 - nPart;
1299 nGenSpecC = 416 - nPartC;
1300 nGenSpecA = 416 - nPartA;
1301 if(nGenSpec>416) nGenSpec=416; if(nGenSpec<0) nGenSpec=0;
1302 if(nGenSpecC>416) nGenSpecC=416; if(nGenSpecC<0) nGenSpecC=0;
1303 if(nGenSpecA>416) nGenSpecA=416; if(nGenSpecA<0) nGenSpecA=0;
1306 delete lineC; delete lineA;
1308 } // ONLY IF fIsCalibrationMB==kFALSE
1310 Bool_t energyFlag = kTRUE;
1311 AliZDCReco* reco = new AliZDCReco(calibSumZN1, calibSumZP1, calibSumZN2, calibSumZP2,
1312 calibTowZN1, calibTowZP1, calibTowZN2, calibTowZP2,
1313 calibZEM1, calibZEM2, sPMRef1, sPMRef2,
1314 nDetSpecNLeft, nDetSpecPLeft, nDetSpecNRight, nDetSpecPRight,
1315 nGenSpec, nGenSpecA, nGenSpecC,
1316 nPart, nPartA, nPartC, b, bA, bC,
1317 recoFlag, energyFlag, isScalerOn, scaler, tdcData);
1319 const Int_t kBufferSize = 4000;
1320 clustersTree->Branch("ZDC", "AliZDCReco", &reco, kBufferSize);
1322 // write the output tree
1323 clustersTree->Fill();
1328 //_____________________________________________________________________________
1329 void AliZDCReconstructor::FillZDCintoESD(TTree *clustersTree, AliESDEvent* esd) const
1331 // fill energies and number of participants to the ESD
1333 // Retrieving TDC calibration data
1334 // Parameters for TDC centering around zero
1335 int const knTDC = 6;
1336 Float_t tdcOffset[knTDC];
1337 for(Int_t jj=0; jj<knTDC; jj++) tdcOffset[jj] = fTDCCalibData->GetMeanTDC(jj);
1338 //fTDCCalibData->Print("");
1341 AliZDCReco* preco = &reco;
1342 clustersTree->SetBranchAddress("ZDC", &preco);
1343 clustersTree->GetEntry(0);
1345 Float_t tZN1Ene[5], tZN2Ene[5], tZP1Ene[5], tZP2Ene[5];
1346 Float_t tZN1EneLR[5], tZN2EneLR[5], tZP1EneLR[5], tZP2EneLR[5];
1347 for(Int_t i=0; i<5; i++){
1348 tZN1Ene[i] = reco.GetZN1HREnTow(i);
1349 tZN2Ene[i] = reco.GetZN2HREnTow(i);
1350 tZP1Ene[i] = reco.GetZP1HREnTow(i);
1351 tZP2Ene[i] = reco.GetZP2HREnTow(i);
1353 tZN1EneLR[i] = reco.GetZN1LREnTow(i);
1354 tZN2EneLR[i] = reco.GetZN2LREnTow(i);
1355 tZP1EneLR[i] = reco.GetZP1LREnTow(i);
1356 tZP2EneLR[i] = reco.GetZP2LREnTow(i);
1359 fESDZDC->SetZN1TowerEnergy(tZN1Ene);
1360 fESDZDC->SetZN2TowerEnergy(tZN2Ene);
1361 fESDZDC->SetZP1TowerEnergy(tZP1Ene);
1362 fESDZDC->SetZP2TowerEnergy(tZP2Ene);
1364 fESDZDC->SetZN1TowerEnergyLR(tZN1EneLR);
1365 fESDZDC->SetZN2TowerEnergyLR(tZN2EneLR);
1366 fESDZDC->SetZP1TowerEnergyLR(tZP1EneLR);
1367 fESDZDC->SetZP2TowerEnergyLR(tZP2EneLR);
1369 Int_t nPart = reco.GetNParticipants();
1370 Int_t nPartA = reco.GetNPartSideA();
1371 Int_t nPartC = reco.GetNPartSideC();
1372 Double_t b = reco.GetImpParameter();
1373 Double_t bA = reco.GetImpParSideA();
1374 Double_t bC = reco.GetImpParSideC();
1375 UInt_t recoFlag = reco.GetRecoFlag();
1377 fESDZDC->SetZDC(reco.GetZN1HREnergy(), reco.GetZP1HREnergy(),
1378 reco.GetZEM1HRsignal(), reco.GetZEM2HRsignal(),
1379 reco.GetZN2HREnergy(), reco.GetZP2HREnergy(),
1380 nPart, nPartA, nPartC, b, bA, bC, recoFlag);
1382 // Writing ZDC scaler for cross section calculation
1383 // ONLY IF the scaler has been read during the event
1384 if(reco.IsScalerOn()==kTRUE){
1386 for(Int_t jk=0; jk<32; jk++) counts[jk] = reco.GetZDCScaler(jk);
1387 fESDZDC->SetZDCScaler(counts);
1390 Int_t tdcValues[32][4] = {{0,}};
1391 Float_t tdcCorrected[32][4] = {{9999.,}};
1392 for(Int_t jk=0; jk<32; jk++){
1393 for(Int_t lk=0; lk<4; lk++){
1394 tdcValues[jk][lk] = reco.GetZDCTDCData(jk, lk);
1396 if(jk==8 && TMath::Abs(tdcValues[jk][lk])>1e-09) fESDZDC->SetZEM1TDChit(kTRUE);
1397 else if(jk==9 && TMath::Abs(tdcValues[jk][lk])>1e-09) fESDZDC->SetZEM2TDChit(kTRUE);
1398 else if(jk==10 && TMath::Abs(tdcValues[jk][lk])>1e-09) fESDZDC->SetZNCTDChit(kTRUE);
1399 else if(jk==11 && TMath::Abs(tdcValues[jk][lk])>1e-09) fESDZDC->SetZPCTDChit(kTRUE);
1400 else if(jk==12 && TMath::Abs(tdcValues[jk][lk])>1e-09) fESDZDC->SetZNATDChit(kTRUE);
1401 else if(jk==13 && TMath::Abs(tdcValues[jk][lk])>1e-09) fESDZDC->SetZPATDChit(kTRUE);
1403 //if((jk>=8 && jk<=13 && lk==0) || jk==15) printf(" *** ZDC: tdc%d = %d = %f ns \n",jk,tdcValues[jk][lk],0.025*tdcValues[jk][lk]);
1407 // Writing TDC data into ZDC ESDs
1408 // 4/2/2011 -> Subtracting L0 (tdcValues[15]) instead of ADC gate
1409 // we try to keep the TDC oscillations as low as possible!
1410 for(Int_t jk=0; jk<32; jk++){
1411 for(Int_t lk=0; lk<4; lk++){
1412 if(tdcValues[jk][lk]!=0.){
1413 // Feb2013 _-> TDC correct entry is there ONLY IF tdc has a hit!
1414 if(TMath::Abs(tdcValues[jk][lk])>1e-09){
1415 tdcCorrected[jk][lk] = 0.025*(tdcValues[jk][lk]-tdcValues[15][0])+fMeanPhase;
1416 // Sep 2011: TDC ch. from 8 to 13 centered around 0 using OCDB
1417 if(jk>=8 && jk<=13) tdcCorrected[jk][lk] = tdcCorrected[jk][lk] - tdcOffset[jk-8];
1419 //if(jk>=8 && jk<=13) printf(" *** tdcOffset%d %f tdcCorr%d %f \n",jk,tdcOffset[jk-8],tdcCorrected[jk][lk]);
1425 fESDZDC->SetZDCTDCData(tdcValues);
1426 fESDZDC->SetZDCTDCCorrected(tdcCorrected);
1427 fESDZDC->AliESDZDC::SetBit(AliESDZDC::kCorrectedTDCFilled, reco.GetEnergyFlag());
1428 fESDZDC->AliESDZDC::SetBit(AliESDZDC::kEnergyCalibratedSignal, kTRUE);
1430 if(esd) esd->SetZDCData(fESDZDC);
1433 //_____________________________________________________________________________
1434 AliCDBStorage* AliZDCReconstructor::SetStorage(const char *uri)
1436 // Setting the storage
1438 Bool_t deleteManager = kFALSE;
1440 AliCDBManager *manager = AliCDBManager::Instance();
1441 AliCDBStorage *defstorage = manager->GetDefaultStorage();
1443 if(!defstorage || !(defstorage->Contains("ZDC"))){
1444 AliWarning("No default storage set or default storage doesn't contain ZDC!");
1445 manager->SetDefaultStorage(uri);
1446 deleteManager = kTRUE;
1449 AliCDBStorage *storage = manager->GetDefaultStorage();
1452 AliCDBManager::Instance()->UnsetDefaultStorage();
1453 defstorage = 0; // the storage is killed by AliCDBManager::Instance()->Destroy()
1459 //_____________________________________________________________________________
1460 AliZDCPedestals* AliZDCReconstructor::GetPedestalData() const
1463 // Getting pedestal calibration object for ZDC set
1465 AliCDBEntry *entry = AliCDBManager::Instance()->Get("ZDC/Calib/Pedestals");
1466 if(!entry) AliFatal("No calibration data loaded!");
1467 entry->SetOwner(kFALSE);
1469 AliZDCPedestals *calibdata = dynamic_cast<AliZDCPedestals*> (entry->GetObject());
1470 if(!calibdata) AliFatal("Wrong calibration object in calibration file!");
1475 //_____________________________________________________________________________
1476 AliZDCEnCalib* AliZDCReconstructor::GetEnergyCalibData() const
1479 // Getting energy and equalization calibration object for ZDC set
1481 AliCDBEntry *entry = AliCDBManager::Instance()->Get("ZDC/Calib/EnergyCalib");
1482 if(!entry) AliFatal("No calibration data loaded!");
1483 entry->SetOwner(kFALSE);
1485 AliZDCEnCalib *calibdata = dynamic_cast<AliZDCEnCalib*> (entry->GetObject());
1486 if(!calibdata) AliFatal("Wrong calibration object in calibration file!");
1491 //_____________________________________________________________________________
1492 AliZDCSaturationCalib* AliZDCReconstructor::GetSaturationCalibData() const
1495 // Getting energy and equalization calibration object for ZDC set
1497 AliCDBEntry *entry = AliCDBManager::Instance()->Get("ZDC/Calib/SaturationCalib");
1498 if(!entry) AliFatal("No calibration data loaded!");
1499 entry->SetOwner(kFALSE);
1501 AliZDCSaturationCalib *calibdata = dynamic_cast<AliZDCSaturationCalib*> (entry->GetObject());
1502 if(!calibdata) AliFatal("Wrong calibration object in calibration file!");
1507 //_____________________________________________________________________________
1508 AliZDCTowerCalib* AliZDCReconstructor::GetTowerCalibData() const
1511 // Getting energy and equalization calibration object for ZDC set
1513 AliCDBEntry *entry = AliCDBManager::Instance()->Get("ZDC/Calib/TowerCalib");
1514 if(!entry) AliFatal("No calibration data loaded!");
1515 entry->SetOwner(kFALSE);
1517 AliZDCTowerCalib *calibdata = dynamic_cast<AliZDCTowerCalib*> (entry->GetObject());
1518 if(!calibdata) AliFatal("Wrong calibration object in calibration file!");
1523 //_____________________________________________________________________________
1524 AliZDCMBCalib* AliZDCReconstructor::GetMBCalibData() const
1527 // Getting energy and equalization calibration object for ZDC set
1529 AliCDBEntry *entry = AliCDBManager::Instance()->Get("ZDC/Calib/MBCalib");
1530 if(!entry) AliFatal("No calibration data loaded!");
1531 entry->SetOwner(kFALSE);
1533 AliZDCMBCalib *calibdata = dynamic_cast<AliZDCMBCalib*> (entry->GetObject());
1534 if(!calibdata) AliFatal("Wrong calibration object in calibration file!");
1539 //_____________________________________________________________________________
1540 AliZDCTDCCalib* AliZDCReconstructor::GetTDCCalibData() const
1543 // Getting TDC object for ZDC
1545 AliCDBEntry *entry = AliCDBManager::Instance()->Get("ZDC/Calib/TDCCalib");
1546 if(!entry) AliFatal("No calibration data loaded!");
1547 entry->SetOwner(kFALSE);
1549 AliZDCTDCCalib *calibdata = dynamic_cast<AliZDCTDCCalib*> (entry->GetObject());
1550 if(!calibdata) AliFatal("Wrong calibration object in calibration file!");