Getting rid of removed class AliPHOSEPFlattener
[u/mrichter/AliRoot.git] / ZDC / AliZDCReconstructor.cxx
CommitLineData
8309c1ab 1/**************************************************************************
2 * Copyright(c) 1998-1999, ALICE Experiment at CERN, All rights reserved. *
3 * *
4 * Author: The ALICE Off-line Project. *
5 * Contributors are mentioned in the code where appropriate. *
6 * *
7 * Permission to use, copy, modify and distribute this software and its *
8 * documentation strictly for non-commercial purposes is hereby granted *
9 * without fee, provided that the above copyright notice appears in all *
10 * copies and that both the copyright notice and this permission notice *
11 * appear in the supporting documentation. The authors make no claims *
12 * about the suitability of this software for any purpose. It is *
f5d41205 13 * provided "as is" without express or implied warranty. *
14 **************************************************************************/
15
16/* $Id$ */
17
18///////////////////////////////////////////////////////////////////////////////
19// //
73bc3a3f 20// ************** Class for ZDC reconstruction ************** //
21// Author: Chiara.Oppedisano@to.infn.it //
22// //
23// NOTATIONS ADOPTED TO IDENTIFY DETECTORS (used in different ages!): //
24// (ZN1,ZP1) or (ZNC, ZPC) or RIGHT refers to side C (RB26) //
25// (ZN2,ZP2) or (ZNA, ZPA) or LEFT refers to side A (RB24) //
f5d41205 26// //
27///////////////////////////////////////////////////////////////////////////////
28
29
73bc3a3f 30#include <TH2F.h>
31#include <TH1D.h>
32#include <TAxis.h>
fd9afd60 33#include <TMap.h>
f5d41205 34
f5d41205 35#include "AliRawReader.h"
36#include "AliESDEvent.h"
a85132e7 37#include "AliESDZDC.h"
f5d41205 38#include "AliZDCDigit.h"
39#include "AliZDCRawStream.h"
40#include "AliZDCReco.h"
41#include "AliZDCReconstructor.h"
6024ec85 42#include "AliZDCPedestals.h"
73bc3a3f 43#include "AliZDCEnCalib.h"
44#include "AliZDCTowerCalib.h"
0d579f58 45#include "AliZDCMBCalib.h"
27d1ff1f 46#include "AliZDCTDCCalib.h"
7bff3766 47#include "AliZDCRecoParam.h"
48#include "AliZDCRecoParampp.h"
49#include "AliZDCRecoParamPbPb.h"
4a72fbdb 50#include "AliRunInfo.h"
9e05925b 51#include "AliLHCClockPhase.h"
f5d41205 52
53
54ClassImp(AliZDCReconstructor)
90936733 55AliZDCRecoParam *AliZDCReconstructor::fgRecoParam=0; //reconstruction parameters
56AliZDCMBCalib *AliZDCReconstructor::fgMBCalibData=0; //calibration parameters for A-A reconstruction
f5d41205 57
58//_____________________________________________________________________________
59AliZDCReconstructor:: AliZDCReconstructor() :
1e319f71 60 fPedData(GetPedestalData()),
61 fEnCalibData(GetEnergyCalibData()),
62 fTowCalibData(GetTowerCalibData()),
27d1ff1f 63 fTDCCalibData(GetTDCCalibData()),
fd9afd60 64 fRecoMode(0),
42d8b8d5 65 fBeamEnergy(0.),
73bc3a3f 66 fNRun(0),
67 fIsCalibrationMB(kFALSE),
68 fPedSubMode(0),
af6f24c9 69 fSignalThreshold(7),
9e05925b 70 fMeanPhase(0),
af6f24c9 71 fESDZDC(NULL)
f5d41205 72{
73 // **** Default constructor
f5d41205 74}
75
76
77//_____________________________________________________________________________
78AliZDCReconstructor::~AliZDCReconstructor()
79{
80// destructor
90936733 81// if(fgRecoParam) delete fgRecoParam;
73bc3a3f 82 if(fPedData) delete fPedData;
83 if(fEnCalibData) delete fEnCalibData;
84 if(fTowCalibData) delete fTowCalibData;
af6f24c9 85 if(fgMBCalibData) delete fgMBCalibData;
86 if(fESDZDC) delete fESDZDC;
f5d41205 87}
88
fd9afd60 89//____________________________________________________________________________
73bc3a3f 90void AliZDCReconstructor::Init()
fd9afd60 91{
9e05925b 92 // Setting reconstruction parameters
93
4a72fbdb 94 TString runType = GetRunInfo()->GetRunType();
4a72fbdb 95 if((runType.CompareTo("CALIBRATION_MB")) == 0){
96 fIsCalibrationMB = kTRUE;
97 }
73bc3a3f 98
4a72fbdb 99 TString beamType = GetRunInfo()->GetBeamType();
100 // This is a temporary solution to allow reconstruction in tests without beam
55d9b744 101 if(((beamType.CompareTo("UNKNOWN"))==0) &&
102 ((runType.CompareTo("PHYSICS"))==0 || (runType.CompareTo("CALIBRATION_BC"))==0)){
4a72fbdb 103 fRecoMode=1;
104 }
81f09162 105 /*else if((beamType.CompareTo("UNKNOWN"))==0){
106 AliError("\t UNKNOWN beam type\n");
4a72fbdb 107 return;
81f09162 108 }*/
9e05925b 109
110 fBeamEnergy = GetRunInfo()->GetBeamEnergy();
dcb6916a 111 if(fBeamEnergy<0.01){
112 AliWarning(" Beam energy value missing -> setting it to 1380 GeV ");
113 fBeamEnergy = 1380.;
114 }
115
4a72fbdb 116 if(((beamType.CompareTo("pp"))==0) || ((beamType.CompareTo("p-p"))==0)
81f09162 117 ||((beamType.CompareTo("PP"))==0) || ((beamType.CompareTo("P-P"))==0)){
118 fRecoMode=1;
119 }
ee4b2f68 120 else if(((beamType.CompareTo("p-A"))==0) || ((beamType.CompareTo("A-p"))==0)
f6c61b40 121 ||((beamType.CompareTo("P-A"))==0) || ((beamType.CompareTo("A-P"))==0)){
122 fRecoMode=1;
123 }
2d9c70ab 124 else if((beamType.CompareTo("A-A")) == 0 || (beamType.CompareTo("AA")) == 0){
4a72fbdb 125 fRecoMode=2;
9e05925b 126 if(!fgRecoParam) fgRecoParam = const_cast<AliZDCRecoParam*>(GetRecoParam());
127 if(fgRecoParam){
128 fgRecoParam->SetGlauberMCDist(fBeamEnergy);
129 }
fd9afd60 130 }
9e05925b 131
132 AliCDBEntry *entry = AliCDBManager::Instance()->Get("GRP/Calib/LHCClockPhase");
133 if (!entry) AliFatal("LHC clock-phase shift is not found in OCDB !");
134 AliLHCClockPhase *phaseLHC = (AliLHCClockPhase*)entry->GetObject();
4162afb1 135 // 4/2/2011 According to A. Di Mauro BEAM1 measurement is more reliable
136 // than BEAM2 and therefore also than the average of the 2
137 fMeanPhase = phaseLHC->GetMeanPhaseB1();
4a72fbdb 138
139 if(fIsCalibrationMB==kFALSE)
b184d3dd 140 AliInfo(Form("\n\n ***** ZDC reconstruction initialized for %s @ %1.0f + %1.0f GeV *****\n\n",
141 beamType.Data(), fBeamEnergy, fBeamEnergy));
fd9afd60 142
32e2fda5 143 // if EMD calibration run NO ENERGY CALIBRATION should be performed
144 // pp-like reconstruction must be performed (E cailb. coeff. = 1)
145 if((runType.CompareTo("CALIBRATION_EMD")) == 0){
146 fRecoMode=1;
9e05925b 147 fBeamEnergy = 1380.;
32e2fda5 148 }
149
b64f2ac3 150 AliInfo(Form("\n ZDC reconstruction mode %d (1 -> p-p, 2-> A-A)\n\n",fRecoMode));
151
af6f24c9 152 fESDZDC = new AliESDZDC();
153
154}
155
156
157//____________________________________________________________________________
158void AliZDCReconstructor::Init(TString beamType, Float_t beamEnergy)
159{
160 // Setting reconstruction mode
161 // Needed to work in the HLT framework
162
163 fIsCalibrationMB = kFALSE;
164
165 fBeamEnergy = beamEnergy;
166
167 if(((beamType.CompareTo("pp"))==0) || ((beamType.CompareTo("p-p"))==0)
168 ||((beamType.CompareTo("PP"))==0) || ((beamType.CompareTo("P-P"))==0)){
169 fRecoMode=1;
170 }
ee4b2f68 171 else if(((beamType.CompareTo("p-A"))==0) || ((beamType.CompareTo("A-p"))==0)
172 ||((beamType.CompareTo("P-A"))==0) || ((beamType.CompareTo("A-P"))==0)){
173 fRecoMode=1;
174 }
af6f24c9 175 else if((beamType.CompareTo("A-A")) == 0 || (beamType.CompareTo("AA")) == 0){
176 fRecoMode=2;
bb98da29 177 if(!fgRecoParam) fgRecoParam = const_cast<AliZDCRecoParam*>(GetRecoParam());
178 if( fgRecoParam ) fgRecoParam->SetGlauberMCDist(fBeamEnergy);
179 }
9e05925b 180
181 AliCDBEntry *entry = AliCDBManager::Instance()->Get("GRP/Calib/LHCClockPhase");
182 if (!entry) AliFatal("LHC clock-phase shift is not found in OCDB !");
183 AliLHCClockPhase *phaseLHC = (AliLHCClockPhase*)entry->GetObject();
184 fMeanPhase = phaseLHC->GetMeanPhase();
af6f24c9 185
186 fESDZDC = new AliESDZDC();
187
b184d3dd 188 AliInfo(Form("\n\n ***** ZDC reconstruction initialized for %s @ %1.0f + %1.0f GeV *****\n\n",
189 beamType.Data(), fBeamEnergy, fBeamEnergy));
af6f24c9 190
fd9afd60 191}
192
f5d41205 193//_____________________________________________________________________________
1e319f71 194void AliZDCReconstructor::Reconstruct(TTree* digitsTree, TTree* clustersTree) const
f5d41205 195{
196 // *** Local ZDC reconstruction for digits
197 // Works on the current event
198
199 // Retrieving calibration data
42d8b8d5 200 // Parameters for mean value pedestal subtraction
73bc3a3f 201 int const kNch = 24;
202 Float_t meanPed[2*kNch];
203 for(Int_t jj=0; jj<2*kNch; jj++) meanPed[jj] = fPedData->GetMeanPed(jj);
42d8b8d5 204 // Parameters pedestal subtraction through correlation with out-of-time signals
73bc3a3f 205 Float_t corrCoeff0[2*kNch], corrCoeff1[2*kNch];
206 for(Int_t jj=0; jj<2*kNch; jj++){
81f09162 207 corrCoeff0[jj] = fPedData->GetPedCorrCoeff0(jj);
208 corrCoeff1[jj] = fPedData->GetPedCorrCoeff1(jj);
42d8b8d5 209 }
f5d41205 210
211 // get digits
212 AliZDCDigit digit;
213 AliZDCDigit* pdigit = &digit;
214 digitsTree->SetBranchAddress("ZDC", &pdigit);
c35ed519 215 //printf("\n\t # of digits in tree: %d\n",(Int_t) digitsTree->GetEntries());
f5d41205 216
217 // loop over digits
c35ed519 218 Float_t tZN1Corr[10], tZP1Corr[10], tZN2Corr[10], tZP2Corr[10];
73bc3a3f 219 Float_t dZEM1Corr[2], dZEM2Corr[2], sPMRef1[2], sPMRef2[2];
c35ed519 220 for(Int_t i=0; i<10; i++){
221 tZN1Corr[i] = tZP1Corr[i] = tZN2Corr[i] = tZP2Corr[i] = 0.;
73bc3a3f 222 if(i<2) dZEM1Corr[i] = dZEM2Corr[i] = sPMRef1[i] = sPMRef2[i] = 0.;
c35ed519 223 }
42d8b8d5 224
225 Int_t digNentries = digitsTree->GetEntries();
796c8b58 226 Float_t ootDigi[kNch]; Int_t i=0;
42d8b8d5 227 // -- Reading out-of-time signals (last kNch entries) for current event
228 if(fPedSubMode==1){
229 for(Int_t iDigit=kNch; iDigit<digNentries; iDigit++){
796c8b58 230 if(i<=kNch) ootDigi[i] = digitsTree->GetEntry(iDigit);
231 else AliWarning(" Can't read more out of time values: index>kNch !!!\n");
232 i++;
42d8b8d5 233 }
234 }
235
236 for(Int_t iDigit=0; iDigit<(digNentries/2); iDigit++) {
a85132e7 237 digitsTree->GetEntry(iDigit);
238 if (!pdigit) continue;
a85132e7 239 //
240 Int_t det = digit.GetSector(0);
241 Int_t quad = digit.GetSector(1);
42d8b8d5 242 Int_t pedindex = -1;
243 Float_t ped2SubHg=0., ped2SubLg=0.;
244 if(quad!=5){
245 if(det==1) pedindex = quad;
246 else if(det==2) pedindex = quad+5;
247 else if(det==3) pedindex = quad+9;
248 else if(det==4) pedindex = quad+12;
249 else if(det==5) pedindex = quad+17;
250 }
251 else pedindex = (det-1)/3+22;
a85132e7 252 //
42d8b8d5 253 if(fPedSubMode==0){
254 ped2SubHg = meanPed[pedindex];
255 ped2SubLg = meanPed[pedindex+kNch];
256 }
257 else if(fPedSubMode==1){
258 ped2SubHg = corrCoeff1[pedindex]*ootDigi[pedindex]+corrCoeff0[pedindex];
259 ped2SubLg = corrCoeff1[pedindex+kNch]*ootDigi[pedindex+kNch]+corrCoeff0[pedindex+kNch];
260 }
42d8b8d5 261
a85132e7 262 if(quad != 5){ // ZDC (not reference PTMs!)
c35ed519 263 if(det == 1){ // *** ZNC
42d8b8d5 264 tZN1Corr[quad] = (Float_t) (digit.GetADCValue(0)-ped2SubHg);
265 tZN1Corr[quad+5] = (Float_t) (digit.GetADCValue(1)-ped2SubLg);
f5d41205 266 }
267 else if(det == 2){ // *** ZP1
42d8b8d5 268 tZP1Corr[quad] = (Float_t) (digit.GetADCValue(0)-ped2SubHg);
42d8b8d5 269 tZP1Corr[quad+5] = (Float_t) (digit.GetADCValue(1)-ped2SubLg);
f5d41205 270 }
271 else if(det == 3){
272 if(quad == 1){ // *** ZEM1
42d8b8d5 273 dZEM1Corr[0] += (Float_t) (digit.GetADCValue(0)-ped2SubHg);
42d8b8d5 274 dZEM1Corr[1] += (Float_t) (digit.GetADCValue(1)-ped2SubLg);
f5d41205 275 }
a85132e7 276 else if(quad == 2){ // *** ZEM2
42d8b8d5 277 dZEM2Corr[0] += (Float_t) (digit.GetADCValue(0)-ped2SubHg);
42d8b8d5 278 dZEM2Corr[1] += (Float_t) (digit.GetADCValue(1)-ped2SubLg);
f5d41205 279 }
280 }
281 else if(det == 4){ // *** ZN2
42d8b8d5 282 tZN2Corr[quad] = (Float_t) (digit.GetADCValue(0)-ped2SubHg);
42d8b8d5 283 tZN2Corr[quad+5] = (Float_t) (digit.GetADCValue(1)-ped2SubLg);
59a953e0 284 }
f5d41205 285 else if(det == 5){ // *** ZP2
42d8b8d5 286 tZP2Corr[quad] = (Float_t) (digit.GetADCValue(0)-ped2SubHg);
42d8b8d5 287 tZP2Corr[quad+5] = (Float_t) (digit.GetADCValue(1)-ped2SubLg);
f5d41205 288 }
a85132e7 289 }
c35ed519 290 else{ // Reference PMs
c35ed519 291 if(det == 1){
73bc3a3f 292 sPMRef1[0] = (Float_t) (digit.GetADCValue(0)-ped2SubHg);
293 sPMRef1[1] = (Float_t) (digit.GetADCValue(1)-ped2SubLg);
c35ed519 294 }
295 else if(det == 4){
73bc3a3f 296 sPMRef2[0] = (Float_t) (digit.GetADCValue(0)-ped2SubHg);
297 sPMRef2[1] = (Float_t) (digit.GetADCValue(1)-ped2SubLg);
c35ed519 298 }
299 }
f5d41205 300
73bc3a3f 301 // Ch. debug
59a953e0 302 /*printf("AliZDCReconstructor: digit #%d det %d quad %d pedHG %1.0f pedLG %1.0f\n",
73bc3a3f 303 iDigit, det, quad, ped2SubHg, ped2SubLg);
59a953e0 304 printf(" -> pedindex %d\n", pedindex);
305 printf(" HGChain -> RawDig %d DigCorr %1.2f",
306 digit.GetADCValue(0), digit.GetADCValue(0)-ped2SubHg);
307 printf(" LGChain -> RawDig %d DigCorr %1.2f\n",
308 digit.GetADCValue(1), digit.GetADCValue(1)-ped2SubLg);*/
309
73bc3a3f 310 }//digits loop
58671297 311
81f09162 312 UInt_t counts[32];
f53e5ecb 313 Int_t tdc[32][4];
82dffa48 314 for(Int_t jj=0; jj<32; jj++){
315 counts[jj]=0;
f53e5ecb 316 for(Int_t ii=0; ii<4; ii++) tdc[jj][ii]=0;
82dffa48 317 }
81f09162 318
319 Int_t evQualityBlock[4] = {1,0,0,0};
320 Int_t triggerBlock[4] = {0,0,0,0};
321 Int_t chBlock[3] = {0,0,0};
322 UInt_t puBits=0;
73bc3a3f 323
69550cf5 324 // reconstruct the event
73bc3a3f 325 if(fRecoMode==1)
fd9afd60 326 ReconstructEventpp(clustersTree, tZN1Corr, tZP1Corr, tZN2Corr, tZP2Corr,
81f09162 327 dZEM1Corr, dZEM2Corr, sPMRef1, sPMRef2,
82dffa48 328 kFALSE, counts, tdc,
81f09162 329 evQualityBlock, triggerBlock, chBlock, puBits);
73bc3a3f 330 else if(fRecoMode==2)
fd9afd60 331 ReconstructEventPbPb(clustersTree, tZN1Corr, tZP1Corr, tZN2Corr, tZP2Corr,
81f09162 332 dZEM1Corr, dZEM2Corr, sPMRef1, sPMRef2,
82dffa48 333 kFALSE, counts, tdc,
91debfd4 334 evQualityBlock, triggerBlock, chBlock, puBits);
f5d41205 335}
336
337//_____________________________________________________________________________
1e319f71 338void AliZDCReconstructor::Reconstruct(AliRawReader* rawReader, TTree* clustersTree) const
f5d41205 339{
340 // *** ZDC raw data reconstruction
341 // Works on the current event
342
343 // Retrieving calibration data
73bc3a3f 344 // Parameters for pedestal subtraction
345 int const kNch = 24;
346 Float_t meanPed[2*kNch];
347 for(Int_t jj=0; jj<2*kNch; jj++) meanPed[jj] = fPedData->GetMeanPed(jj);
42d8b8d5 348 // Parameters pedestal subtraction through correlation with out-of-time signals
73bc3a3f 349 Float_t corrCoeff0[2*kNch], corrCoeff1[2*kNch];
350 for(Int_t jj=0; jj<2*kNch; jj++){
42d8b8d5 351 corrCoeff0[jj] = fPedData->GetPedCorrCoeff0(jj);
352 corrCoeff1[jj] = fPedData->GetPedCorrCoeff1(jj);
36660aad 353 //printf(" %d %1.4f %1.4f\n", jj,corrCoeff0[jj],corrCoeff1[jj]);
42d8b8d5 354 }
f5d41205 355
73bc3a3f 356 Int_t adcZN1[5], adcZN1oot[5], adcZN1lg[5], adcZN1ootlg[5];
357 Int_t adcZP1[5], adcZP1oot[5], adcZP1lg[5], adcZP1ootlg[5];
358 Int_t adcZN2[5], adcZN2oot[5], adcZN2lg[5], adcZN2ootlg[5];
359 Int_t adcZP2[5], adcZP2oot[5], adcZP2lg[5], adcZP2ootlg[5];
360 Int_t adcZEM[2], adcZEMoot[2], adcZEMlg[2], adcZEMootlg[2];
361 Int_t pmRef[2], pmRefoot[2], pmReflg[2], pmRefootlg[2];
362 for(Int_t ich=0; ich<5; ich++){
363 adcZN1[ich] = adcZN1oot[ich] = adcZN1lg[ich] = adcZN1ootlg[ich] = 0;
364 adcZP1[ich] = adcZP1oot[ich] = adcZP1lg[ich] = adcZP1ootlg[ich] = 0;
365 adcZN2[ich] = adcZN2oot[ich] = adcZN2lg[ich] = adcZN2ootlg[ich] = 0;
366 adcZP2[ich] = adcZP2oot[ich] = adcZP2lg[ich] = adcZP2ootlg[ich] = 0;
367 if(ich<2){
368 adcZEM[ich] = adcZEMoot[ich] = adcZEMlg[ich] = adcZEMootlg[ich] = 0;
369 pmRef[ich] = pmRefoot[ich] = pmReflg[ich] = pmRefootlg[ich] = 0;
370 }
371 }
7bff3766 372
c35ed519 373 Float_t tZN1Corr[10], tZP1Corr[10], tZN2Corr[10], tZP2Corr[10];
73bc3a3f 374 Float_t dZEM1Corr[2], dZEM2Corr[2], sPMRef1[2], sPMRef2[2];
c35ed519 375 for(Int_t i=0; i<10; i++){
376 tZN1Corr[i] = tZP1Corr[i] = tZN2Corr[i] = tZP2Corr[i] = 0.;
73bc3a3f 377 if(i<2) dZEM1Corr[i] = dZEM2Corr[i] = sPMRef1[i] = sPMRef2[i] = 0.;
c35ed519 378 }
1e319f71 379
81f09162 380 Bool_t isScalerOn=kFALSE;
32e2fda5 381 Int_t jsc=0, itdc=0, iprevtdc=-1, ihittdc=0;
82dffa48 382 UInt_t scalerData[32];
f53e5ecb 383 Int_t tdcData[32][4];
82dffa48 384 for(Int_t k=0; k<32; k++){
385 scalerData[k]=0;
010f62f2 386 for(Int_t i=0; i<4; i++) tdcData[k][i]=0;
82dffa48 387 }
81f09162 388
9e05925b 389
81f09162 390 Int_t evQualityBlock[4] = {1,0,0,0};
391 Int_t triggerBlock[4] = {0,0,0,0};
392 Int_t chBlock[3] = {0,0,0};
393 UInt_t puBits=0;
394
82dffa48 395 Int_t kFirstADCGeo=0, kLastADCGeo=3, kScalerGeo=8, kZDCTDCGeo=4, kPUGeo=29;
81f09162 396 //Int_t kTrigScales=30, kTrigHistory=31;
397
398 // loop over raw data
af6f24c9 399 //rawReader->Reset();
f5d41205 400 AliZDCRawStream rawData(rawReader);
73bc3a3f 401 while(rawData.Next()){
81f09162 402
403 // ***************************** Reading ADCs
f70a5526 404 if((rawData.GetADCModule()>=kFirstADCGeo) && (rawData.GetADCModule()<=kLastADCGeo)){
405 //printf(" **** Reading ADC raw data from module %d **** \n",rawData.GetADCModule());
58671297 406 //
81f09162 407 if((rawData.IsADCDataWord()) && (rawData.GetNChannelsOn()<48)) chBlock[0] = kTRUE;
408 if((rawData.IsADCDataWord()) && (rawData.IsOverflow() == kTRUE)) chBlock[1] = kTRUE;
409 if((rawData.IsADCDataWord()) && (rawData.IsUnderflow() == kTRUE)) chBlock[2] = kTRUE;
410 if((rawData.IsADCDataWord()) && (rawData.IsADCEventGood() == kTRUE)) evQualityBlock[0] = kTRUE;
f70a5526 411
81f09162 412 if((rawData.IsADCDataWord()) && (rawData.IsUnderflow()==kFALSE)
413 && (rawData.IsOverflow()==kFALSE) && (rawData.IsADCEventGood()==kTRUE)){
73bc3a3f 414
58671297 415 Int_t adcMod = rawData.GetADCModule();
416 Int_t det = rawData.GetSector(0);
417 Int_t quad = rawData.GetSector(1);
418 Int_t gain = rawData.GetADCGain();
419 Int_t pedindex=0;
420 //
421 // Mean pedestal value subtraction -------------------------------------------------------
422 if(fPedSubMode == 0){
36660aad 423 // **** Pb-Pb data taking 2010 -> subtracting some ch. from correlation ****
73bc3a3f 424 // Not interested in o.o.t. signals (ADC modules 2, 3)
36660aad 425 //if(adcMod == 2 || adcMod == 3) continue;
b64f2ac3 426 // **** Pb-Pb data taking 2011 -> subtracting only ZEM from correlation ****
427 if(det==3){
428 if(adcMod==0 || adcMod==1){
429 if(gain==0) adcZEM[quad-1] = rawData.GetADCValue();
430 else adcZEMlg[quad-1] = rawData.GetADCValue();
36660aad 431 }
b64f2ac3 432 else if(adcMod==2 || adcMod==3){
433 if(gain==0) adcZEMoot[quad-1] = rawData.GetADCValue();
434 else adcZEMootlg[quad-1] = rawData.GetADCValue();
36660aad 435 }
436 }
437 // When oot values are read the ADC modules 2, 3 can be skipped!!!
1e319f71 438 if(adcMod == 2 || adcMod == 3) continue;
36660aad 439
440 // *************************************************************************
73bc3a3f 441 if(quad != 5){ // ZDCs (not reference PTMs)
b64f2ac3 442 if(det==1){
73bc3a3f 443 pedindex = quad;
444 if(gain == 0) tZN1Corr[quad] += (Float_t) (rawData.GetADCValue()-meanPed[pedindex]);
445 else tZN1Corr[quad+5] += (Float_t) (rawData.GetADCValue()-meanPed[pedindex+kNch]);
f5d41205 446 }
077d2505 447 else if(det==2){
73bc3a3f 448 pedindex = quad+5;
449 if(gain == 0) tZP1Corr[quad] += (Float_t) (rawData.GetADCValue()-meanPed[pedindex]);
450 else tZP1Corr[quad+5] += (Float_t) (rawData.GetADCValue()-meanPed[pedindex+kNch]);
451 }
36660aad 452 /*else if(det == 3){
73bc3a3f 453 pedindex = quad+9;
454 if(quad==1){
455 if(gain == 0) dZEM1Corr[0] += (Float_t) (rawData.GetADCValue()-meanPed[pedindex]);
456 else dZEM1Corr[1] += (Float_t) (rawData.GetADCValue()-meanPed[pedindex+kNch]);
457 }
458 else if(quad==2){
459 if(gain == 0) dZEM2Corr[0] += (Float_t) (rawData.GetADCValue()-meanPed[pedindex]);
460 else dZEM2Corr[1] += (Float_t) (rawData.GetADCValue()-meanPed[pedindex+kNch]);
461 }
36660aad 462 }*/
73bc3a3f 463 else if(det == 4){
464 pedindex = quad+12;
465 if(gain == 0) tZN2Corr[quad] += (Float_t) (rawData.GetADCValue()-meanPed[pedindex]);
466 else tZN2Corr[quad+5] += (Float_t) (rawData.GetADCValue()-meanPed[pedindex+kNch]);
467 }
468 else if(det == 5){
469 pedindex = quad+17;
470 if(gain == 0) tZP2Corr[quad] += (Float_t) (rawData.GetADCValue()-meanPed[pedindex]);
471 else tZP2Corr[quad+5] += (Float_t) (rawData.GetADCValue()-meanPed[pedindex+kNch]);
f5d41205 472 }
c35ed519 473 }
73bc3a3f 474 else{ // reference PM
475 pedindex = (det-1)/3 + 22;
476 if(det == 1){
477 if(gain==0) sPMRef1[0] += (Float_t) (rawData.GetADCValue()-meanPed[pedindex]);
81f09162 478 else sPMRef1[1] += (Float_t) (rawData.GetADCValue()-meanPed[pedindex+kNch]);
73bc3a3f 479 }
480 else if(det == 4){
481 if(gain==0) sPMRef2[0] += (Float_t) (rawData.GetADCValue()-meanPed[pedindex]);
81f09162 482 else sPMRef2[1] += (Float_t) (rawData.GetADCValue()-meanPed[pedindex+kNch]);
73bc3a3f 483 }
c35ed519 484 }
73bc3a3f 485 // Ch. debug
81f09162 486 /*if(gain==0){
487 printf(" AliZDCReconstructor: det %d quad %d res %d -> Pedestal[%d] %1.0f",
488 det,quad,gain, pedindex, meanPed[pedindex]);
489 printf(" RawADC %d ADCCorr %1.0f\n",
490 rawData.GetADCValue(), rawData.GetADCValue()-meanPed[pedindex]);
491 }*/
58671297 492 }// mean pedestal subtraction
493 // Pedestal subtraction from correlation ------------------------------------------------
494 else if(fPedSubMode == 1){
73bc3a3f 495 // In time signals
496 if(adcMod==0 || adcMod==1){
497 if(quad != 5){ // signals from ZDCs
498 if(det == 1){
499 if(gain==0) adcZN1[quad] = rawData.GetADCValue();
500 else adcZN1lg[quad] = rawData.GetADCValue();
501 }
502 else if(det == 2){
503 if(gain==0) adcZP1[quad] = rawData.GetADCValue();
504 else adcZP1lg[quad] = rawData.GetADCValue();
505 }
506 else if(det == 3){
507 if(gain==0) adcZEM[quad-1] = rawData.GetADCValue();
508 else adcZEMlg[quad-1] = rawData.GetADCValue();
509 }
510 else if(det == 4){
511 if(gain==0) adcZN2[quad] = rawData.GetADCValue();
512 else adcZN2lg[quad] = rawData.GetADCValue();
513 }
514 else if(det == 5){
515 if(gain==0) adcZP2[quad] = rawData.GetADCValue();
516 else adcZP2lg[quad] = rawData.GetADCValue();
517 }
518 }
519 else{ // signals from reference PM
520 if(gain==0) pmRef[quad-1] = rawData.GetADCValue();
521 else pmReflg[quad-1] = rawData.GetADCValue();
522 }
523 }
524 // Out-of-time pedestals
525 else if(adcMod==2 || adcMod==3){
526 if(quad != 5){ // signals from ZDCs
527 if(det == 1){
528 if(gain==0) adcZN1oot[quad] = rawData.GetADCValue();
529 else adcZN1ootlg[quad] = rawData.GetADCValue();
530 }
531 else if(det == 2){
532 if(gain==0) adcZP1oot[quad] = rawData.GetADCValue();
533 else adcZP1ootlg[quad] = rawData.GetADCValue();
534 }
535 else if(det == 3){
536 if(gain==0) adcZEMoot[quad-1] = rawData.GetADCValue();
537 else adcZEMootlg[quad-1] = rawData.GetADCValue();
538 }
539 else if(det == 4){
540 if(gain==0) adcZN2oot[quad] = rawData.GetADCValue();
541 else adcZN2ootlg[quad] = rawData.GetADCValue();
542 }
543 else if(det == 5){
544 if(gain==0) adcZP2oot[quad] = rawData.GetADCValue();
545 else adcZP2ootlg[quad] = rawData.GetADCValue();
546 }
547 }
548 else{ // signals from reference PM
549 if(gain==0) pmRefoot[quad-1] = rawData.GetADCValue();
550 else pmRefootlg[quad-1] = rawData.GetADCValue();
551 }
552 }
58671297 553 } // pedestal subtraction from correlation
554 // Ch. debug
59a953e0 555 /*printf("\t AliZDCReconstructor: det %d quad %d res %d -> Ped[%d] = %1.0f\n",
556 det,quad,gain, pedindex, meanPed[pedindex]);*/
58671297 557 }//IsADCDataWord
81f09162 558 }// ADC DATA
559 // ***************************** Reading Scaler
560 else if(rawData.GetADCModule()==kScalerGeo){
27d1ff1f 561 if(rawData.IsScalerWord()==kTRUE){
81f09162 562 isScalerOn = kTRUE;
563 scalerData[jsc] = rawData.GetTriggerCount();
ad3a602e 564 // Ch. debug
565 //printf(" Reconstructed VME Scaler: %d %d ",jsc,scalerData[jsc]);
566 //
81f09162 567 jsc++;
568 }
569 }// VME SCALER DATA
82dffa48 570 // ***************************** Reading ZDC TDC
571 else if(rawData.GetADCModule()==kZDCTDCGeo && rawData.IsZDCTDCDatum()==kTRUE){
32e2fda5 572 itdc = rawData.GetChannel();
573 if(itdc==iprevtdc) ihittdc++;
574 else ihittdc=0;
575 iprevtdc=itdc;
650731d2 576 if(ihittdc<4) tdcData[itdc][ihittdc] = rawData.GetZDCTDCDatum();
82dffa48 577 // Ch. debug
27d1ff1f 578 //if(ihittdc==0) printf(" TDC%d %d ",itdc, tdcData[itdc][ihittdc]);
82dffa48 579 }// ZDC TDC DATA
81f09162 580 // ***************************** Reading PU
581 else if(rawData.GetADCModule()==kPUGeo){
582 puBits = rawData.GetDetectorPattern();
58671297 583 }
81f09162 584 // ***************************** Reading trigger history
585 else if(rawData.IstriggerHistoryWord()==kTRUE){
586 triggerBlock[0] = rawData.IsCPTInputEMDTrigger();
587 triggerBlock[1] = rawData.IsCPTInputSemiCentralTrigger();
588 triggerBlock[2] = rawData.IsCPTInputCentralTrigger();
589 triggerBlock[3] = rawData.IsCPTInputMBTrigger();
590 }
591
73bc3a3f 592 }//loop on raw data
593
594 if(fPedSubMode==1){
595 for(Int_t t=0; t<5; t++){
596 tZN1Corr[t] = adcZN1[t] - (corrCoeff1[t]*adcZN1oot[t]+corrCoeff0[t]);
597 tZN1Corr[t+5] = adcZN1lg[t] - (corrCoeff1[t+kNch]*adcZN1ootlg[t]+corrCoeff0[t+kNch]);
598 //
599 tZP1Corr[t] = adcZP1[t] - (corrCoeff1[t+5]*adcZP1oot[t]+corrCoeff0[t+5]);
600 tZP1Corr[t+5] = adcZP1lg[t] - (corrCoeff1[t+5+kNch]*adcZP1ootlg[t]+corrCoeff0[t+5+kNch]);
601 //
602 tZN2Corr[t] = adcZN2[t] - (corrCoeff1[t+12]*adcZN2oot[t]+corrCoeff0[t+12]);
603 tZN2Corr[t+5] = adcZN2lg[t] - (corrCoeff1[t+12+kNch]*adcZN2ootlg[t]+corrCoeff0[t+12+kNch]);
604 //
605 tZP2Corr[t] = adcZP2[t] - (corrCoeff1[t+17]*adcZP2oot[t]+corrCoeff0[t+17]);
606 tZP2Corr[t+5] = adcZP2lg[t] - (corrCoeff1[t+17+kNch]*adcZP2ootlg[t]+corrCoeff0[t+17+kNch]);
73bc3a3f 607 }
36660aad 608 dZEM1Corr[0] = adcZEM[0] - (corrCoeff1[10]*adcZEMoot[0]+corrCoeff0[10]);
609 dZEM1Corr[1] = adcZEMlg[0] - (corrCoeff1[10+kNch]*adcZEMootlg[0]+corrCoeff0[10+kNch]);
610 dZEM2Corr[0] = adcZEM[1] - (corrCoeff1[11]*adcZEMoot[1]+corrCoeff0[11]);
611 dZEM2Corr[1] = adcZEMlg[1] - (corrCoeff1[11+kNch]*adcZEMootlg[1]+corrCoeff0[11+kNch]);
73bc3a3f 612 //
613 sPMRef1[0] = pmRef[0] - (corrCoeff1[22]*pmRefoot[0]+corrCoeff0[22]);
614 sPMRef1[1] = pmReflg[0] - (corrCoeff1[22+kNch]*pmRefootlg[0]+corrCoeff0[22+kNch]);
615 sPMRef2[0] = pmRef[0] - (corrCoeff1[23]*pmRefoot[1]+corrCoeff0[23]);
616 sPMRef2[1] = pmReflg[0] - (corrCoeff1[23+kNch]*pmRefootlg[1]+corrCoeff0[23+kNch]);
f5d41205 617 }
ae7346df 618 if(fPedSubMode==0 && fRecoMode==2){
b64f2ac3 619 // **** Pb-Pb data taking 2011 -> subtracting some ch. from correlation ****
620 //tZN1Corr[0] = adcZN1[0] - (corrCoeff1[0]*adcZN1oot[0]+corrCoeff0[0]);
621 //tZN1Corr[5] = adcZN1lg[0] - (corrCoeff1[kNch]*adcZN1ootlg[0]+corrCoeff0[kNch]);
36660aad 622 // Ch. debug
623 //printf(" adcZN1 %d adcZN1oot %d tZN1Corr %1.2f \n", adcZN1[0],adcZN1oot[0],tZN1Corr[0]);
624 //printf(" adcZN1lg %d adcZN1ootlg %d tZN1Corrlg %1.2f \n", adcZN1lg[0],adcZN1ootlg[0],tZN1Corr[5]);
625 //
077d2505 626 //tZP1Corr[2] = adcZP1[2] - (corrCoeff1[2+5]*adcZP1oot[2]+corrCoeff0[2+5]);
627 //tZP1Corr[2+5] = adcZP1lg[2] - (corrCoeff1[2+5+kNch]*adcZP1ootlg[2]+corrCoeff0[2+5+kNch]);
36660aad 628 //
629 dZEM1Corr[0] = adcZEM[0] - (corrCoeff1[10]*adcZEMoot[0]+corrCoeff0[10]);
630 dZEM1Corr[1] = adcZEMlg[0] - (corrCoeff1[10+kNch]*adcZEMootlg[0]+corrCoeff0[10+kNch]);
631 dZEM2Corr[0] = adcZEM[1] - (corrCoeff1[11]*adcZEMoot[1]+corrCoeff0[11]);
632 dZEM2Corr[1] = adcZEMlg[1] - (corrCoeff1[11+kNch]*adcZEMootlg[1]+corrCoeff0[11+kNch]);
633 // *************************************************************************
634 }
ae7346df 635 else if(fPedSubMode==0 && fRecoMode==1){
636 // **** p-p data taking 2011 -> temporary patch to overcome DA problem ****
ae7346df 637 //
638 dZEM1Corr[0] = adcZEM[0] - meanPed[10];
639 dZEM1Corr[1] = adcZEMlg[0] - meanPed[10+kNch];
640 dZEM2Corr[0] = adcZEM[1] - meanPed[11];
641 dZEM2Corr[1] = adcZEMlg[1] - meanPed[11+kNch];
642 // *************************************************************************
643 }
f5d41205 644
81f09162 645 if(fRecoMode==1) // p-p data
646 ReconstructEventpp(clustersTree, tZN1Corr, tZP1Corr, tZN2Corr, tZP2Corr,
647 dZEM1Corr, dZEM2Corr, sPMRef1, sPMRef2,
82dffa48 648 isScalerOn, scalerData, tdcData,
81f09162 649 evQualityBlock, triggerBlock, chBlock, puBits);
650 else if(fRecoMode==2) // Pb-Pb data
91debfd4 651 ReconstructEventPbPb(clustersTree, tZN1Corr, tZP1Corr, tZN2Corr, tZP2Corr,
81f09162 652 dZEM1Corr, dZEM2Corr, sPMRef1, sPMRef2,
82dffa48 653 isScalerOn, scalerData, tdcData,
81f09162 654 evQualityBlock, triggerBlock, chBlock, puBits);
f5d41205 655}
656
657//_____________________________________________________________________________
90936733 658void AliZDCReconstructor::ReconstructEventpp(TTree *clustersTree,
659 const Float_t* const corrADCZN1, const Float_t* const corrADCZP1,
660 const Float_t* const corrADCZN2, const Float_t* const corrADCZP2,
661 const Float_t* const corrADCZEM1, const Float_t* const corrADCZEM2,
662 Float_t* sPMRef1, Float_t* sPMRef2, Bool_t isScalerOn, UInt_t* scaler,
f53e5ecb 663 Int_t tdcData[32][4], const Int_t* const evQualityBlock,
82dffa48 664 const Int_t* const triggerBlock, const Int_t* const chBlock, UInt_t puBits) const
f5d41205 665{
73bc3a3f 666 // ****************** Reconstruct one event ******************
59a953e0 667
668 // CH. debug
669 /*printf("\n*************************************************\n");
670 printf(" ReconstructEventpp -> values after pedestal subtraction:\n");
671 printf(" ADCZN1 [%1.2f %1.2f %1.2f %1.2f %1.2f]\n",
672 corrADCZN1[0],corrADCZN1[1],corrADCZN1[2],corrADCZN1[3],corrADCZN1[4]);
673 printf(" ADCZP1 [%1.2f %1.2f %1.2f %1.2f %1.2f]\n",
674 corrADCZP1[0],corrADCZP1[1],corrADCZP1[2],corrADCZP1[3],corrADCZP1[4]);
675 printf(" ADCZN2 [%1.2f %1.2f %1.2f %1.2f %1.2f]\n",
676 corrADCZN2[0],corrADCZN2[1],corrADCZN2[2],corrADCZN2[3],corrADCZN2[4]);
677 printf(" ADCZP2 [%1.2f %1.2f %1.2f %1.2f %1.2f]\n",
678 corrADCZP2[0],corrADCZP2[1],corrADCZP2[2],corrADCZP2[3],corrADCZP2[4]);
679 printf(" ADCZEM1 [%1.2f] ADCZEM2 [%1.2f] \n",corrADCZEM1[0],corrADCZEM2[0]);
680 printf("*************************************************\n");*/
1e319f71 681
81f09162 682 // ---------------------- Setting reco flags for ESD
6b793021 683 UInt_t rFlags[32];
684 for(Int_t ifl=0; ifl<32; ifl++) rFlags[ifl]=0;
81f09162 685
686 if(evQualityBlock[0] == 1) rFlags[31] = 0x0;
687 else rFlags[31] = 0x1;
1e319f71 688 //
81f09162 689 if(evQualityBlock[1] == 1) rFlags[30] = 0x1;
690 if(evQualityBlock[2] == 1) rFlags[29] = 0x1;
691 if(evQualityBlock[3] == 1) rFlags[28] = 0x1;
692
693 if(triggerBlock[0] == 1) rFlags[27] = 0x1;
694 if(triggerBlock[1] == 1) rFlags[26] = 0x1;
695 if(triggerBlock[2] == 1) rFlags[25] = 0x1;
696 if(triggerBlock[3] == 1) rFlags[24] = 0x1;
697
698 if(chBlock[0] == 1) rFlags[18] = 0x1;
699 if(chBlock[1] == 1) rFlags[17] = 0x1;
700 if(chBlock[2] == 1) rFlags[16] = 0x1;
701
702
703 rFlags[13] = puBits & 0x00000020;
704 rFlags[12] = puBits & 0x00000010;
705 rFlags[11] = puBits & 0x00000080;
706 rFlags[10] = puBits & 0x00000040;
707 rFlags[9] = puBits & 0x00000020;
708 rFlags[8] = puBits & 0x00000010;
709
710 if(corrADCZP1[0]>fSignalThreshold) rFlags[5] = 0x1;
711 if(corrADCZN1[0]>fSignalThreshold) rFlags[4] = 0x1;
712 if(corrADCZEM2[0]>fSignalThreshold) rFlags[3] = 0x1;
713 if(corrADCZEM1[0]>fSignalThreshold) rFlags[2] = 0x1;
714 if(corrADCZP2[0]>fSignalThreshold) rFlags[1] = 0x1;
715 if(corrADCZN2[0]>fSignalThreshold) rFlags[0] = 0x1;
716
717 UInt_t recoFlag = rFlags[31] << 31 | rFlags[30] << 30 | rFlags[29] << 29 | rFlags[28] << 28 |
718 rFlags[27] << 27 | rFlags[26] << 26 | rFlags[25] << 25 | rFlags[24] << 24 |
719 0x0 << 23 | 0x0 << 22 | 0x0 << 21 | 0x0 << 20 |
720 0x0 << 19 | rFlags[18] << 18 | rFlags[17] << 17 | rFlags[16] << 16 |
721 0x0 << 15 | 0x0 << 14 | rFlags[13] << 13 | rFlags[12] << 12 |
722 rFlags[11] << 11 |rFlags[10] << 10 | rFlags[9] << 9 | rFlags[8] << 8 |
723 0x0 << 7 | 0x0 << 6 | rFlags[5] << 5 | rFlags[4] << 4 |
724 rFlags[3] << 3 | rFlags[2] << 2 | rFlags[1] << 1 | rFlags[0];
725 // --------------------------------------------------
1e319f71 726
73bc3a3f 727 // ****** Retrieving calibration data
84d6255e 728 // --- Equalization coefficients ---------------------------------------------
f5d41205 729 Float_t equalCoeffZN1[5], equalCoeffZP1[5], equalCoeffZN2[5], equalCoeffZP2[5];
730 for(Int_t ji=0; ji<5; ji++){
73bc3a3f 731 equalCoeffZN1[ji] = fTowCalibData->GetZN1EqualCoeff(ji);
732 equalCoeffZP1[ji] = fTowCalibData->GetZP1EqualCoeff(ji);
733 equalCoeffZN2[ji] = fTowCalibData->GetZN2EqualCoeff(ji);
734 equalCoeffZP2[ji] = fTowCalibData->GetZP2EqualCoeff(ji);
f5d41205 735 }
84d6255e 736 // --- Energy calibration factors ------------------------------------
796c8b58 737 Float_t calibEne[6];
42d8b8d5 738 // **** Energy calibration coefficient set to 1
739 // **** (no trivial way to calibrate in p-p runs)
76725070 740 for(Int_t ij=0; ij<6; ij++) calibEne[ij] = fEnCalibData->GetEnCalib(ij);
7bff3766 741
73bc3a3f 742 // ****** Equalization of detector responses
7bff3766 743 Float_t equalTowZN1[10], equalTowZN2[10], equalTowZP1[10], equalTowZP2[10];
73bc3a3f 744 for(Int_t gi=0; gi<10; gi++){
31474197 745 if(gi<5){
746 equalTowZN1[gi] = corrADCZN1[gi]*equalCoeffZN1[gi];
747 equalTowZP1[gi] = corrADCZP1[gi]*equalCoeffZP1[gi];
748 equalTowZN2[gi] = corrADCZN2[gi]*equalCoeffZN2[gi];
749 equalTowZP2[gi] = corrADCZP2[gi]*equalCoeffZP2[gi];
750 }
751 else{
752 equalTowZN1[gi] = corrADCZN1[gi]*equalCoeffZN1[gi-5];
753 equalTowZP1[gi] = corrADCZP1[gi]*equalCoeffZP1[gi-5];
754 equalTowZN2[gi] = corrADCZN2[gi]*equalCoeffZN2[gi-5];
755 equalTowZP2[gi] = corrADCZP2[gi]*equalCoeffZP2[gi-5];
756 }
73bc3a3f 757 }
59a953e0 758 // Ch. debug
759 /*printf("\n ------------- EQUALIZATION -------------\n");
760 printf(" ADCZN1 [%1.2f %1.2f %1.2f %1.2f %1.2f]\n",
761 equalTowZN1[0],equalTowZN1[1],equalTowZN1[2],equalTowZN1[3],equalTowZN1[4]);
762 printf(" ADCZP1 [%1.2f %1.2f %1.2f %1.2f %1.2f]\n",
763 equalTowZP1[0],equalTowZP1[1],equalTowZP1[2],equalTowZP1[3],equalTowZP1[4]);
764 printf(" ADCZN2 [%1.2f %1.2f %1.2f %1.2f %1.2f]\n",
765 equalTowZN2[0],equalTowZN2[1],equalTowZN2[2],equalTowZN2[3],equalTowZN2[4]);
766 printf(" ADCZP2 [%1.2f %1.2f %1.2f %1.2f %1.2f]\n",
767 equalTowZP2[0],equalTowZP2[1],equalTowZP2[2],equalTowZP2[3],equalTowZP2[4]);
768 printf(" ----------------------------------------\n");*/
73bc3a3f 769
770 // ****** Summed response for hadronic calorimeter (SUMMED and then CALIBRATED!)
771 Float_t calibSumZN1[]={0,0}, calibSumZN2[]={0,0}, calibSumZP1[]={0,0}, calibSumZP2[]={0,0};
7bff3766 772 for(Int_t gi=0; gi<5; gi++){
73bc3a3f 773 calibSumZN1[0] += equalTowZN1[gi];
774 calibSumZP1[0] += equalTowZP1[gi];
775 calibSumZN2[0] += equalTowZN2[gi];
776 calibSumZP2[0] += equalTowZP2[gi];
777 //
778 calibSumZN1[1] += equalTowZN1[gi+5];
779 calibSumZP1[1] += equalTowZP1[gi+5];
780 calibSumZN2[1] += equalTowZN2[gi+5];
781 calibSumZP2[1] += equalTowZP2[gi+5];
7bff3766 782 }
73bc3a3f 783 // High gain chain
59a953e0 784 calibSumZN1[0] = calibSumZN1[0]*calibEne[0];
785 calibSumZP1[0] = calibSumZP1[0]*calibEne[1];
786 calibSumZN2[0] = calibSumZN2[0]*calibEne[2];
787 calibSumZP2[0] = calibSumZP2[0]*calibEne[3];
73bc3a3f 788 // Low gain chain
789 calibSumZN1[1] = calibSumZN1[1]*calibEne[0];
790 calibSumZP1[1] = calibSumZP1[1]*calibEne[1];
791 calibSumZN2[1] = calibSumZN2[1]*calibEne[2];
792 calibSumZP2[1] = calibSumZP2[1]*calibEne[3];
7bff3766 793
73bc3a3f 794 // ****** Energy calibration of detector responses
7bff3766 795 Float_t calibTowZN1[10], calibTowZN2[10], calibTowZP1[10], calibTowZP2[10];
73bc3a3f 796 for(Int_t gi=0; gi<5; gi++){
797 // High gain chain
59a953e0 798 calibTowZN1[gi] = equalTowZN1[gi]*calibEne[0];
799 calibTowZP1[gi] = equalTowZP1[gi]*calibEne[1];
800 calibTowZN2[gi] = equalTowZN2[gi]*calibEne[2];
801 calibTowZP2[gi] = equalTowZP2[gi]*calibEne[3];
73bc3a3f 802 // Low gain chain
803 calibTowZN1[gi+5] = equalTowZN1[gi+5]*calibEne[0];
804 calibTowZP1[gi+5] = equalTowZP1[gi+5]*calibEne[1];
805 calibTowZN2[gi+5] = equalTowZN2[gi+5]*calibEne[2];
806 calibTowZP2[gi+5] = equalTowZP2[gi+5]*calibEne[3];
7bff3766 807 }
73bc3a3f 808 //
809 Float_t sumZEM[]={0,0}, calibZEM1[]={0,0}, calibZEM2[]={0,0};
c7d86465 810 calibZEM1[0] = corrADCZEM1[0]*calibEne[4];
811 calibZEM1[1] = corrADCZEM1[1]*calibEne[4];
59a953e0 812 calibZEM2[0] = corrADCZEM2[0]*calibEne[5];
73bc3a3f 813 calibZEM2[1] = corrADCZEM2[1]*calibEne[5];
814 for(Int_t k=0; k<2; k++) sumZEM[k] = calibZEM1[k] + calibZEM2[k];
59a953e0 815 // Ch. debug
816 /*printf("\n ------------- CALIBRATION -------------\n");
817 printf(" ADCZN1 [%1.2f %1.2f %1.2f %1.2f %1.2f]\n",
818 calibTowZN1[0],calibTowZN1[1],calibTowZN1[2],calibTowZN1[3],calibTowZN1[4]);
819 printf(" ADCZP1 [%1.2f %1.2f %1.2f %1.2f %1.2f]\n",
820 calibTowZP1[0],calibTowZP1[1],calibTowZP1[2],calibTowZP1[3],calibTowZP1[4]);
821 printf(" ADCZN2 [%1.2f %1.2f %1.2f %1.2f %1.2f]\n",
822 calibTowZN2[0],calibTowZN2[1],calibTowZN2[2],calibTowZN2[3],calibTowZN2[4]);
823 printf(" ADCZP2 [%1.2f %1.2f %1.2f %1.2f %1.2f]\n",
824 calibTowZP2[0],calibTowZP2[1],calibTowZP2[2],calibTowZP2[3],calibTowZP2[4]);
825 printf(" ADCZEM1 [%1.2f] ADCZEM2 [%1.2f] \n",calibZEM1[0],calibZEM2[0]);
826 printf(" ----------------------------------------\n");*/
7bff3766 827
73bc3a3f 828 // ****** No. of spectator and participants nucleons
42d8b8d5 829 // Variables calculated to comply with ESD structure
73bc3a3f 830 // *** N.B. -> They have a meaning only in Pb-Pb!!!!!!!!!!!!
d9ec113e 831 Int_t nDetSpecNLeft=0, nDetSpecPLeft=0, nDetSpecNRight=0, nDetSpecPRight=0;
73bc3a3f 832 Int_t nGenSpec=0, nGenSpecLeft=0, nGenSpecRight=0;
833 Int_t nPart=0, nPartTotLeft=0, nPartTotRight=0;
834 Double_t impPar=0., impPar1=0., impPar2=0.;
7bff3766 835
077d2505 836 Bool_t energyFlag = kFALSE;
7bff3766 837 // create the output tree
1e319f71 838 AliZDCReco* reco = new AliZDCReco(calibSumZN1, calibSumZP1, calibSumZN2, calibSumZP2,
839 calibTowZN1, calibTowZP1, calibTowZN2, calibTowZP2,
840 calibZEM1, calibZEM2, sPMRef1, sPMRef2,
841 nDetSpecNLeft, nDetSpecPLeft, nDetSpecNRight, nDetSpecPRight,
842 nGenSpec, nGenSpecLeft, nGenSpecRight,
843 nPart, nPartTotLeft, nPartTotRight,
844 impPar, impPar1, impPar2,
077d2505 845 recoFlag, energyFlag, isScalerOn, scaler, tdcData);
7bff3766 846
1b467dea 847 const Int_t kBufferSize = 4000;
1e319f71 848 clustersTree->Branch("ZDC", "AliZDCReco", &reco, kBufferSize);
7bff3766 849 // write the output tree
850 clustersTree->Fill();
bb98da29 851 delete reco;
7bff3766 852}
853
854//_____________________________________________________________________________
73bc3a3f 855void AliZDCReconstructor::ReconstructEventPbPb(TTree *clustersTree,
90936733 856 const Float_t* const corrADCZN1, const Float_t* const corrADCZP1,
857 const Float_t* const corrADCZN2, const Float_t* const corrADCZP2,
858 const Float_t* const corrADCZEM1, const Float_t* const corrADCZEM2,
859 Float_t* sPMRef1, Float_t* sPMRef2, Bool_t isScalerOn, UInt_t* scaler,
f53e5ecb 860 Int_t tdcData[32][4], const Int_t* const evQualityBlock,
82dffa48 861 const Int_t* const triggerBlock, const Int_t* const chBlock, UInt_t puBits) const
7bff3766 862{
73bc3a3f 863 // ****************** Reconstruct one event ******************
81f09162 864 // ---------------------- Setting reco flags for ESD
6b793021 865 UInt_t rFlags[32];
866 for(Int_t ifl=0; ifl<32; ifl++) rFlags[ifl]=0;
81f09162 867
868 if(evQualityBlock[0] == 1) rFlags[31] = 0x0;
869 else rFlags[31] = 0x1;
1e319f71 870 //
81f09162 871 if(evQualityBlock[1] == 1) rFlags[30] = 0x1;
872 if(evQualityBlock[2] == 1) rFlags[29] = 0x1;
873 if(evQualityBlock[3] == 1) rFlags[28] = 0x1;
874
875 if(triggerBlock[0] == 1) rFlags[27] = 0x1;
876 if(triggerBlock[1] == 1) rFlags[26] = 0x1;
877 if(triggerBlock[2] == 1) rFlags[25] = 0x1;
878 if(triggerBlock[3] == 1) rFlags[24] = 0x1;
879
880 if(chBlock[0] == 1) rFlags[18] = 0x1;
881 if(chBlock[1] == 1) rFlags[17] = 0x1;
882 if(chBlock[2] == 1) rFlags[16] = 0x1;
883
884 rFlags[13] = puBits & 0x00000020;
885 rFlags[12] = puBits & 0x00000010;
886 rFlags[11] = puBits & 0x00000080;
887 rFlags[10] = puBits & 0x00000040;
888 rFlags[9] = puBits & 0x00000020;
889 rFlags[8] = puBits & 0x00000010;
890
891 if(corrADCZP1[0]>fSignalThreshold) rFlags[5] = 0x1;
892 if(corrADCZN1[0]>fSignalThreshold) rFlags[4] = 0x1;
893 if(corrADCZEM2[0]>fSignalThreshold) rFlags[3] = 0x1;
894 if(corrADCZEM1[0]>fSignalThreshold) rFlags[2] = 0x1;
895 if(corrADCZP2[0]>fSignalThreshold) rFlags[1] = 0x1;
896 if(corrADCZN2[0]>fSignalThreshold) rFlags[0] = 0x1;
897
898 UInt_t recoFlag = rFlags[31] << 31 | rFlags[30] << 30 | rFlags[29] << 29 | rFlags[28] << 28 |
899 rFlags[27] << 27 | rFlags[26] << 26 | rFlags[25] << 25 | rFlags[24] << 24 |
900 0x0 << 23 | 0x0 << 22 | 0x0 << 21 | 0x0 << 20 |
901 0x0 << 19 | rFlags[18] << 18 | rFlags[17] << 17 | rFlags[16] << 16 |
902 0x0 << 15 | 0x0 << 14 | rFlags[13] << 13 | rFlags[12] << 12 |
903 rFlags[11] << 11 |rFlags[10] << 10 | rFlags[9] << 9 | rFlags[8] << 8 |
904 0x0 << 7 | 0x0 << 6 | rFlags[5] << 5 | rFlags[4] << 4 |
905 rFlags[3] << 3 | rFlags[2] << 2 | rFlags[1] << 1 | rFlags[0];
906 // --------------------------------------------------
907
2c62f191 908
909 // CH. debug
910/* printf("\n*************************************************\n");
911 printf(" ReconstructEventPbPb -> values after pedestal subtraction:\n");
912 printf(" ADCZN1 [%1.2f %1.2f %1.2f %1.2f %1.2f]\n",
913 corrADCZN1[0],corrADCZN1[1],corrADCZN1[2],corrADCZN1[3],corrADCZN1[4]);
914 printf(" ADCZP1 [%1.2f %1.2f %1.2f %1.2f %1.2f]\n",
915 corrADCZP1[0],corrADCZP1[1],corrADCZP1[2],corrADCZP1[3],corrADCZP1[4]);
916 printf(" ADCZN2 [%1.2f %1.2f %1.2f %1.2f %1.2f]\n",
917 corrADCZN2[0],corrADCZN2[1],corrADCZN2[2],corrADCZN2[3],corrADCZN2[4]);
918 printf(" ADCZP2 [%1.2f %1.2f %1.2f %1.2f %1.2f]\n",
919 corrADCZP2[0],corrADCZP2[1],corrADCZP2[2],corrADCZP2[3],corrADCZP2[4]);
920 printf(" ADCZEM1 [%1.2f] ADCZEM2 [%1.2f] \n",corrADCZEM1[0],corrADCZEM2[0]);
921 printf("*************************************************\n");
922*/
73bc3a3f 923 // ****** Retrieving calibration data
7bff3766 924 // --- Equalization coefficients ---------------------------------------------
925 Float_t equalCoeffZN1[5], equalCoeffZP1[5], equalCoeffZN2[5], equalCoeffZP2[5];
926 for(Int_t ji=0; ji<5; ji++){
73bc3a3f 927 equalCoeffZN1[ji] = fTowCalibData->GetZN1EqualCoeff(ji);
928 equalCoeffZP1[ji] = fTowCalibData->GetZP1EqualCoeff(ji);
929 equalCoeffZN2[ji] = fTowCalibData->GetZN2EqualCoeff(ji);
930 equalCoeffZP2[ji] = fTowCalibData->GetZP2EqualCoeff(ji);
7bff3766 931 }
932 // --- Energy calibration factors ------------------------------------
2c62f191 933 Float_t calibEne[6];
934 // The energy calibration object already takes into account of E_beam
935 // -> the value from the OCDB can be directly used (Jul 2010)
936 for(Int_t ij=0; ij<6; ij++) calibEne[ij] = fEnCalibData->GetEnCalib(ij);
f5d41205 937
73bc3a3f 938 // ****** Equalization of detector responses
c35ed519 939 Float_t equalTowZN1[10], equalTowZN2[10], equalTowZP1[10], equalTowZP2[10];
73bc3a3f 940 for(Int_t gi=0; gi<10; gi++){
31474197 941 if(gi<5){
942 equalTowZN1[gi] = corrADCZN1[gi]*equalCoeffZN1[gi];
943 equalTowZP1[gi] = corrADCZP1[gi]*equalCoeffZP1[gi];
944 equalTowZN2[gi] = corrADCZN2[gi]*equalCoeffZN2[gi];
945 equalTowZP2[gi] = corrADCZP2[gi]*equalCoeffZP2[gi];
946 }
947 else{
948 equalTowZN1[gi] = corrADCZN1[gi]*equalCoeffZN1[gi-5];
949 equalTowZP1[gi] = corrADCZP1[gi]*equalCoeffZP1[gi-5];
950 equalTowZN2[gi] = corrADCZN2[gi]*equalCoeffZN2[gi-5];
951 equalTowZP2[gi] = corrADCZP2[gi]*equalCoeffZP2[gi-5];
952 }
f5d41205 953 }
954
2c62f191 955 // Ch. debug
956/* printf("\n ------------- EQUALIZATION -------------\n");
957 printf(" ADCZN1 [%1.2f %1.2f %1.2f %1.2f %1.2f]\n",
958 equalTowZN1[0],equalTowZN1[1],equalTowZN1[2],equalTowZN1[3],equalTowZN1[4]);
959 printf(" ADCZP1 [%1.2f %1.2f %1.2f %1.2f %1.2f]\n",
960 equalTowZP1[0],equalTowZP1[1],equalTowZP1[2],equalTowZP1[3],equalTowZP1[4]);
961 printf(" ADCZN2 [%1.2f %1.2f %1.2f %1.2f %1.2f]\n",
962 equalTowZN2[0],equalTowZN2[1],equalTowZN2[2],equalTowZN2[3],equalTowZN2[4]);
963 printf(" ADCZP2 [%1.2f %1.2f %1.2f %1.2f %1.2f]\n",
964 equalTowZP2[0],equalTowZP2[1],equalTowZP2[2],equalTowZP2[3],equalTowZP2[4]);
965 printf(" ----------------------------------------\n");
966*/
967
73bc3a3f 968 // ****** Summed response for hadronic calorimeter (SUMMED and then CALIBRATED!)
c9102a72 969 Float_t calibSumZN1[]={0,0}, calibSumZN2[]={0,0}, calibSumZP1[]={0,0}, calibSumZP2[]={0,0};
73bc3a3f 970 for(Int_t gi=0; gi<5; gi++){
971 calibSumZN1[0] += equalTowZN1[gi];
972 calibSumZP1[0] += equalTowZP1[gi];
973 calibSumZN2[0] += equalTowZN2[gi];
974 calibSumZP2[0] += equalTowZP2[gi];
975 //
976 calibSumZN1[1] += equalTowZN1[gi+5];
977 calibSumZP1[1] += equalTowZP1[gi+5];
978 calibSumZN2[1] += equalTowZN2[gi+5];
979 calibSumZP2[1] += equalTowZP2[gi+5];
f5d41205 980 }
2c62f191 981 //
2413fd18 982 //fEnCalibData->Print("");
2c62f191 983
73bc3a3f 984 // High gain chain
1b32dcc6 985 calibSumZN1[0] = calibSumZN1[0]*calibEne[0]*8.;
986 calibSumZP1[0] = calibSumZP1[0]*calibEne[1]*8.;
987 calibSumZN2[0] = calibSumZN2[0]*calibEne[2]*8.;
988 calibSumZP2[0] = calibSumZP2[0]*calibEne[3]*8.;
73bc3a3f 989 // Low gain chain
990 calibSumZN1[1] = calibSumZN1[1]*calibEne[0];
991 calibSumZP1[1] = calibSumZP1[1]*calibEne[1];
992 calibSumZN2[1] = calibSumZN2[1]*calibEne[2];
993 calibSumZP2[1] = calibSumZP2[1]*calibEne[3];
7bff3766 994 //
73bc3a3f 995 Float_t sumZEM[]={0,0}, calibZEM1[]={0,0}, calibZEM2[]={0,0};
c7d86465 996 calibZEM1[0] = corrADCZEM1[0]*calibEne[4]*8.;
997 calibZEM1[1] = corrADCZEM1[1]*calibEne[4];
1b32dcc6 998 calibZEM2[0] = corrADCZEM2[0]*calibEne[5]*8.;
73bc3a3f 999 calibZEM2[1] = corrADCZEM2[1]*calibEne[5];
1000 for(Int_t k=0; k<2; k++) sumZEM[k] = calibZEM1[k] + calibZEM2[k];
1001
1002 // ****** Energy calibration of detector responses
1003 Float_t calibTowZN1[10], calibTowZN2[10], calibTowZP1[10], calibTowZP2[10];
1004 for(Int_t gi=0; gi<5; gi++){
1005 // High gain chain
6f915a8f 1006 calibTowZN1[gi] = equalTowZN1[gi]*2*calibEne[0]*8.;
1007 calibTowZP1[gi] = equalTowZP1[gi]*2*calibEne[1]*8.;
1008 calibTowZN2[gi] = equalTowZN2[gi]*2*calibEne[2]*8.;
1009 calibTowZP2[gi] = equalTowZP2[gi]*2*calibEne[3]*8.;
73bc3a3f 1010 // Low gain chain
6f915a8f 1011 calibTowZN1[gi+5] = equalTowZN1[gi+5]*2*calibEne[0];
1012 calibTowZP1[gi+5] = equalTowZP1[gi+5]*2*calibEne[1];
1013 calibTowZN2[gi+5] = equalTowZN2[gi+5]*2*calibEne[2];
1014 calibTowZP2[gi+5] = equalTowZP2[gi+5]*2*calibEne[3];
73bc3a3f 1015 }
2c62f191 1016
1017 // Ch. debug
1018/* printf("\n ------------- CALIBRATION -------------\n");
1019 printf(" ADCZN1 [%1.2f %1.2f %1.2f %1.2f %1.2f]\n",
1020 calibTowZN1[0],calibTowZN1[1],calibTowZN1[2],calibTowZN1[3],calibTowZN1[4]);
1021 printf(" ADCZP1 [%1.2f %1.2f %1.2f %1.2f %1.2f]\n",
1022 calibTowZP1[0],calibTowZP1[1],calibTowZP1[2],calibTowZP1[3],calibTowZP1[4]);
1023 printf(" ADCZN2 [%1.2f %1.2f %1.2f %1.2f %1.2f]\n",
1024 calibTowZN2[0],calibTowZN2[1],calibTowZN2[2],calibTowZN2[3],calibTowZN2[4]);
1025 printf(" ADCZP2 [%1.2f %1.2f %1.2f %1.2f %1.2f]\n",
1026 calibTowZP2[0],calibTowZP2[1],calibTowZP2[2],calibTowZP2[3],calibTowZP2[4]);
1027 printf(" ADCZEM1 [%1.2f] ADCZEM2 [%1.2f] \n",calibZEM1[0],calibZEM2[0]);
1028 printf(" ----------------------------------------\n");
1029*/
73bc3a3f 1030 // ****** Number of detected spectator nucleons
d9ec113e 1031 Int_t nDetSpecNLeft=0, nDetSpecPLeft=0, nDetSpecNRight=0, nDetSpecPRight=0;
90936733 1032 if(fBeamEnergy>0.01){
fd9afd60 1033 nDetSpecNLeft = (Int_t) (calibSumZN1[0]/fBeamEnergy);
1034 nDetSpecPLeft = (Int_t) (calibSumZP1[0]/fBeamEnergy);
1035 nDetSpecNRight = (Int_t) (calibSumZN2[0]/fBeamEnergy);
1036 nDetSpecPRight = (Int_t) (calibSumZP2[0]/fBeamEnergy);
1037 }
73bc3a3f 1038 else AliWarning(" ATTENTION!!! fBeamEnergy=0 -> N_spec will be ZERO!!! \n");
af6f24c9 1039 /*printf("\n\t AliZDCReconstructor -> fBeamEnergy %1.0f: nDetSpecNsideA %d, nDetSpecPsideA %d,"
2c62f191 1040 " nDetSpecNsideC %d, nDetSpecPsideC %d\n",fBeamEnergy,nDetSpecNLeft, nDetSpecPLeft,
af6f24c9 1041 nDetSpecNRight, nDetSpecPRight);*/
73bc3a3f 1042
1e319f71 1043 Int_t nGenSpec=0, nGenSpecA=0, nGenSpecC=0;
1044 Int_t nPart=0, nPartA=0, nPartC=0;
1045 Double_t b=0., bA=0., bC=0.;
1046
73bc3a3f 1047 if(fIsCalibrationMB == kFALSE){
6318fa5a 1048 // ****** Reconstruction parameters ------------------
1049 if(!fgRecoParam) fgRecoParam = const_cast<AliZDCRecoParam*>(GetRecoParam());
1050 if(!fgRecoParam){
1051 AliError(" RecoParam object not retrieved correctly: not reconstructing ZDC event!!!");
1052 return;
1053 }
1054 TH1D* hNpartDist = fgRecoParam->GethNpartDist();
1055 TH1D* hbDist = fgRecoParam->GethbDist();
1056 Float_t fClkCenter = fgRecoParam->GetClkCenter();
1057 if(!hNpartDist || !hbDist){
1058 AliError("Something wrong in Glauber MC histos got from AliZDCREcoParamPbPb: NO EVENT RECO FOR ZDC DATA!!!\n\n");
1059 //return;
1060 }
1061 else{
9e05925b 1062 if(!fgMBCalibData) fgMBCalibData = const_cast<AliZDCMBCalib*>(GetMBCalibData());
90936733 1063 TH2F *hZDCvsZEM = fgMBCalibData->GethZDCvsZEM();
1064 TH2F *hZDCCvsZEM = fgMBCalibData->GethZDCCvsZEM();
1065 TH2F *hZDCAvsZEM = fgMBCalibData->GethZDCAvsZEM();
0d579f58 1066 //
73bc3a3f 1067 Double_t xHighEdge = hZDCvsZEM->GetXaxis()->GetXmax();
9e05925b 1068 Double_t origin = xHighEdge*fClkCenter;
73bc3a3f 1069 // Ch. debug
1e319f71 1070 //printf("\n\n xHighEdge %1.2f, origin %1.4f \n", xHighEdge, origin);
73bc3a3f 1071 //
1072 // ====> Summed ZDC info (sideA+side C)
1073 TF1 *line = new TF1("line","[0]*x+[1]",0.,xHighEdge);
1074 Float_t y = (calibSumZN1[0]+calibSumZP1[0]+calibSumZN2[0]+calibSumZP2[0])/1000.;
1075 Float_t x = (calibZEM1[0]+calibZEM2[0])/1000.;
1076 line->SetParameter(0, y/(x-origin));
1077 line->SetParameter(1, -origin*y/(x-origin));
1078 // Ch. debug
1e319f71 1079 //printf(" ***************** Summed ZDC info (sideA+side C) \n");
1080 //printf(" E_{ZEM} %1.4f, E_{ZDC} %1.2f, TF1: %1.2f*x + %1.2f ", x, y,y/(x-origin),-origin*y/(x-origin));
73bc3a3f 1081 //
1082 Double_t countPerc=0;
1083 Double_t xBinCenter=0, yBinCenter=0;
1084 for(Int_t nbinx=1; nbinx<=hZDCvsZEM->GetNbinsX(); nbinx++){
1085 for(Int_t nbiny=1; nbiny<=hZDCvsZEM->GetNbinsY(); nbiny++){
1086 xBinCenter = hZDCvsZEM->GetXaxis()->GetBinCenter(nbinx);
1087 yBinCenter = hZDCvsZEM->GetYaxis()->GetBinCenter(nbiny);
1088 //
1089 if(line->GetParameter(0)>0){
1090 if(yBinCenter < (line->GetParameter(0)*xBinCenter + line->GetParameter(1))){
1091 countPerc += hZDCvsZEM->GetBinContent(nbinx,nbiny);
1092 // Ch. debug
2c62f191 1093 //printf(" xBinCenter %1.3f, yBinCenter %1.0f, countPerc %1.0f\n",
1094 //xBinCenter, yBinCenter, countPerc);
73bc3a3f 1095 }
1096 }
1097 else{
1098 if(yBinCenter > (line->GetParameter(0)*xBinCenter + line->GetParameter(1))){
1099 countPerc += hZDCvsZEM->GetBinContent(nbinx,nbiny);
1100 // Ch. debug
2c62f191 1101 //printf(" xBinCenter %1.3f, yBinCenter %1.0f, countPerc %1.0f\n",
1102 //xBinCenter, yBinCenter, countPerc);
73bc3a3f 1103 }
1104 }
1105 }
1106 }
1107 //
1108 Double_t xSecPerc = 0.;
1109 if(hZDCvsZEM->GetEntries()!=0){
1110 xSecPerc = countPerc/hZDCvsZEM->GetEntries();
1111 }
1112 else{
1113 AliWarning(" Histogram hZDCvsZEM from OCDB has no entries!!!");
1114 }
1115 // Ch. debug
1116 //printf(" xSecPerc %1.4f \n", xSecPerc);
f5d41205 1117
73bc3a3f 1118 // ====> side C
1119 TF1 *lineC = new TF1("lineC","[0]*x+[1]",0.,xHighEdge);
1120 Float_t yC = (calibSumZN1[0]+calibSumZP1[0])/1000.;
1121 lineC->SetParameter(0, yC/(x-origin));
1122 lineC->SetParameter(1, -origin*yC/(x-origin));
1123 // Ch. debug
1124 //printf(" ***************** Side C \n");
1125 //printf(" E_{ZEM} %1.4f, E_{ZDCC} %1.2f, TF1: %1.2f*x + %1.2f ", x, yC,yC/(x-origin),-origin*yC/(x-origin));
1126 //
1127 Double_t countPercC=0;
1128 Double_t xBinCenterC=0, yBinCenterC=0;
1129 for(Int_t nbinx=1; nbinx<=hZDCCvsZEM->GetNbinsX(); nbinx++){
1130 for(Int_t nbiny=1; nbiny<=hZDCCvsZEM->GetNbinsY(); nbiny++){
1131 xBinCenterC = hZDCCvsZEM->GetXaxis()->GetBinCenter(nbinx);
1132 yBinCenterC = hZDCCvsZEM->GetYaxis()->GetBinCenter(nbiny);
1133 if(lineC->GetParameter(0)>0){
1134 if(yBinCenterC < (lineC->GetParameter(0)*xBinCenterC + lineC->GetParameter(1))){
1135 countPercC += hZDCCvsZEM->GetBinContent(nbinx,nbiny);
1136 }
1137 }
1138 else{
1139 if(yBinCenterC > (lineC->GetParameter(0)*xBinCenterC + lineC->GetParameter(1))){
1140 countPercC += hZDCCvsZEM->GetBinContent(nbinx,nbiny);
1141 }
1142 }
1143 }
1144 }
1145 //
1146 Double_t xSecPercC = 0.;
1147 if(hZDCCvsZEM->GetEntries()!=0){
1148 xSecPercC = countPercC/hZDCCvsZEM->GetEntries();
1149 }
1150 else{
1151 AliWarning(" Histogram hZDCCvsZEM from OCDB has no entries!!!");
1152 }
1153 // Ch. debug
1154 //printf(" xSecPercC %1.4f \n", xSecPercC);
1155
1156 // ====> side A
1157 TF1 *lineA = new TF1("lineA","[0]*x+[1]",0.,xHighEdge);
1158 Float_t yA = (calibSumZN2[0]+calibSumZP2[0])/1000.;
1159 lineA->SetParameter(0, yA/(x-origin));
1160 lineA->SetParameter(1, -origin*yA/(x-origin));
1161 //
1162 // Ch. debug
1163 //printf(" ***************** Side A \n");
1164 //printf(" E_{ZEM} %1.4f, E_{ZDCA} %1.2f, TF1: %1.2f*x + %1.2f ", x, yA,yA/(x-origin),-origin*yA/(x-origin));
1165 //
1166 Double_t countPercA=0;
1167 Double_t xBinCenterA=0, yBinCenterA=0;
1168 for(Int_t nbinx=1; nbinx<=hZDCAvsZEM->GetNbinsX(); nbinx++){
1169 for(Int_t nbiny=1; nbiny<=hZDCAvsZEM->GetNbinsY(); nbiny++){
1170 xBinCenterA = hZDCAvsZEM->GetXaxis()->GetBinCenter(nbinx);
1171 yBinCenterA = hZDCAvsZEM->GetYaxis()->GetBinCenter(nbiny);
1172 if(lineA->GetParameter(0)>0){
1173 if(yBinCenterA < (lineA->GetParameter(0)*xBinCenterA + lineA->GetParameter(1))){
1174 countPercA += hZDCAvsZEM->GetBinContent(nbinx,nbiny);
1175 }
1176 }
1177 else{
1178 if(yBinCenterA > (lineA->GetParameter(0)*xBinCenterA + lineA->GetParameter(1))){
1179 countPercA += hZDCAvsZEM->GetBinContent(nbinx,nbiny);
1180 }
1181 }
1182 }
1183 }
1184 //
1185 Double_t xSecPercA = 0.;
1186 if(hZDCAvsZEM->GetEntries()!=0){
1187 xSecPercA = countPercA/hZDCAvsZEM->GetEntries();
1188 }
1189 else{
1190 AliWarning(" Histogram hZDCAvsZEM from OCDB has no entries!!!");
1191 }
1192 // Ch. debug
1193 //printf(" xSecPercA %1.4f \n", xSecPercA);
1194
1195 // ****** Number of participants (from E_ZDC vs. E_ZEM correlation)
73bc3a3f 1196 Double_t nPartFrac=0., nPartFracC=0., nPartFracA=0.;
1197 for(Int_t npbin=1; npbin<hNpartDist->GetNbinsX(); npbin++){
1198 nPartFrac += (hNpartDist->GetBinContent(npbin))/(hNpartDist->GetEntries());
1199 if((1.-nPartFrac) < xSecPerc){
1200 nPart = (Int_t) hNpartDist->GetBinLowEdge(npbin);
1201 // Ch. debug
1202 //printf(" ***************** Summed ZDC info (sideA+side C) \n");
1203 //printf(" nPartFrac %1.4f, nPart %d\n", nPartFrac, nPart);
1204 break;
1205 }
1206 }
1207 if(nPart<0) nPart=0;
1208 //
1209 for(Int_t npbin=1; npbin<hNpartDist->GetNbinsX(); npbin++){
1210 nPartFracC += (hNpartDist->GetBinContent(npbin))/(hNpartDist->GetEntries());
1211 if((1.-nPartFracC) < xSecPercC){
1212 nPartC = (Int_t) hNpartDist->GetBinLowEdge(npbin);
1213 // Ch. debug
1214 //printf(" ***************** Side C \n");
1215 //printf(" nPartFracC %1.4f, nPartC %d\n", nPartFracC, nPartC);
1216 break;
1217 }
1218 }
1219 if(nPartC<0) nPartC=0;
1220 //
1221 for(Int_t npbin=1; npbin<hNpartDist->GetNbinsX(); npbin++){
1222 nPartFracA += (hNpartDist->GetBinContent(npbin))/(hNpartDist->GetEntries());
1223 if((1.-nPartFracA) < xSecPercA){
1224 nPartA = (Int_t) hNpartDist->GetBinLowEdge(npbin);
1225 // Ch. debug
1226 //printf(" ***************** Side A \n");
1227 //printf(" nPartFracA %1.4f, nPartA %d\n\n", nPartFracA, nPartA);
1228 break;
1229 }
1230 }
1231 if(nPartA<0) nPartA=0;
1232
1233 // ****** Impact parameter (from E_ZDC vs. E_ZEM correlation)
73bc3a3f 1234 Double_t bFrac=0., bFracC=0., bFracA=0.;
1235 for(Int_t ibbin=1; ibbin<hbDist->GetNbinsX(); ibbin++){
1236 bFrac += (hbDist->GetBinContent(ibbin))/(hbDist->GetEntries());
e2ec52da 1237 if(bFrac > xSecPerc){
73bc3a3f 1238 b = hbDist->GetBinLowEdge(ibbin);
1239 break;
1240 }
1241 }
1242 //
1243 for(Int_t ibbin=1; ibbin<hbDist->GetNbinsX(); ibbin++){
1244 bFracC += (hbDist->GetBinContent(ibbin))/(hbDist->GetEntries());
e2ec52da 1245 if(bFracC > xSecPercC){
73bc3a3f 1246 bC = hbDist->GetBinLowEdge(ibbin);
1247 break;
1248 }
1249 }
1250 //
1251 for(Int_t ibbin=1; ibbin<hbDist->GetNbinsX(); ibbin++){
e2ec52da 1252 bFracA += (hbDist->GetBinContent(ibbin))/(hbDist->GetEntries());
1253 if(bFracA > xSecPercA){
73bc3a3f 1254 bA = hbDist->GetBinLowEdge(ibbin);
1255 break;
1256 }
1257 }
1258
1259 // ****** Number of spectator nucleons
73bc3a3f 1260 nGenSpec = 416 - nPart;
1261 nGenSpecC = 416 - nPartC;
1262 nGenSpecA = 416 - nPartA;
1263 if(nGenSpec>416) nGenSpec=416; if(nGenSpec<0) nGenSpec=0;
1264 if(nGenSpecC>416) nGenSpecC=416; if(nGenSpecC<0) nGenSpecC=0;
2c62f191 1265 if(nGenSpecA>416) nGenSpecA=416; if(nGenSpecA<0) nGenSpecA=0;
bb98da29 1266
1267 delete line;
73bc3a3f 1268 delete lineC; delete lineA;
6318fa5a 1269 }
73bc3a3f 1270 } // ONLY IF fIsCalibrationMB==kFALSE
1e319f71 1271
077d2505 1272 Bool_t energyFlag = kTRUE;
1e319f71 1273 AliZDCReco* reco = new AliZDCReco(calibSumZN1, calibSumZP1, calibSumZN2, calibSumZP2,
1274 calibTowZN1, calibTowZP1, calibTowZN2, calibTowZP2,
1275 calibZEM1, calibZEM2, sPMRef1, sPMRef2,
1276 nDetSpecNLeft, nDetSpecPLeft, nDetSpecNRight, nDetSpecPRight,
1277 nGenSpec, nGenSpecA, nGenSpecC,
1278 nPart, nPartA, nPartC, b, bA, bC,
077d2505 1279 recoFlag, energyFlag, isScalerOn, scaler, tdcData);
1e319f71 1280
1281 const Int_t kBufferSize = 4000;
1282 clustersTree->Branch("ZDC", "AliZDCReco", &reco, kBufferSize);
c2bb665a 1283 //reco->Print("");
1e319f71 1284 // write the output tree
1285 clustersTree->Fill();
bb98da29 1286 delete reco;
73bc3a3f 1287}
8309c1ab 1288
8309c1ab 1289
1290//_____________________________________________________________________________
70f04f6d 1291void AliZDCReconstructor::FillZDCintoESD(TTree *clustersTree, AliESDEvent* esd) const
8309c1ab 1292{
70f04f6d 1293 // fill energies and number of participants to the ESD
8309c1ab 1294
27d1ff1f 1295 // Retrieving TDC calibration data
1296 // Parameters for TDC centering around zero
1297 int const knTDC = 6;
1298 Float_t tdcOffset[knTDC];
1299 for(Int_t jj=0; jj<knTDC; jj++) tdcOffset[jj] = fTDCCalibData->GetMeanTDC(jj);
1300 //fTDCCalibData->Print("");
1301
8309c1ab 1302 AliZDCReco reco;
1303 AliZDCReco* preco = &reco;
70f04f6d 1304 clustersTree->SetBranchAddress("ZDC", &preco);
70f04f6d 1305 clustersTree->GetEntry(0);
84d6255e 1306 //
a85132e7 1307 Float_t tZN1Ene[5], tZN2Ene[5], tZP1Ene[5], tZP2Ene[5];
1308 Float_t tZN1EneLR[5], tZN2EneLR[5], tZP1EneLR[5], tZP2EneLR[5];
1309 for(Int_t i=0; i<5; i++){
c35ed519 1310 tZN1Ene[i] = reco.GetZN1HREnTow(i);
1311 tZN2Ene[i] = reco.GetZN2HREnTow(i);
1312 tZP1Ene[i] = reco.GetZP1HREnTow(i);
1313 tZP2Ene[i] = reco.GetZP2HREnTow(i);
1314 //
1315 tZN1EneLR[i] = reco.GetZN1LREnTow(i);
1316 tZN2EneLR[i] = reco.GetZN2LREnTow(i);
1317 tZP1EneLR[i] = reco.GetZP1LREnTow(i);
1318 tZP2EneLR[i] = reco.GetZP2LREnTow(i);
e90a5fef 1319 }
73bc3a3f 1320 //
af6f24c9 1321 fESDZDC->SetZN1TowerEnergy(tZN1Ene);
1322 fESDZDC->SetZN2TowerEnergy(tZN2Ene);
1323 fESDZDC->SetZP1TowerEnergy(tZP1Ene);
1324 fESDZDC->SetZP2TowerEnergy(tZP2Ene);
73bc3a3f 1325 //
af6f24c9 1326 fESDZDC->SetZN1TowerEnergyLR(tZN1EneLR);
1327 fESDZDC->SetZN2TowerEnergyLR(tZN2EneLR);
1328 fESDZDC->SetZP1TowerEnergyLR(tZP1EneLR);
1329 fESDZDC->SetZP2TowerEnergyLR(tZP2EneLR);
a85132e7 1330 //
73bc3a3f 1331 Int_t nPart = reco.GetNParticipants();
1332 Int_t nPartA = reco.GetNPartSideA();
1333 Int_t nPartC = reco.GetNPartSideC();
1334 Double_t b = reco.GetImpParameter();
1335 Double_t bA = reco.GetImpParSideA();
1336 Double_t bC = reco.GetImpParSideC();
1e319f71 1337 UInt_t recoFlag = reco.GetRecoFlag();
af6f24c9 1338
1339 fESDZDC->SetZDC(reco.GetZN1HREnergy(), reco.GetZP1HREnergy(),
1340 reco.GetZEM1HRsignal(), reco.GetZEM2HRsignal(),
1341 reco.GetZN2HREnergy(), reco.GetZP2HREnergy(),
1342 nPart, nPartA, nPartC, b, bA, bC, recoFlag);
a4cab348 1343
81f09162 1344 // Writing ZDC scaler for cross section calculation
1345 // ONLY IF the scaler has been read during the event
1346 if(reco.IsScalerOn()==kTRUE){
1347 UInt_t counts[32];
1348 for(Int_t jk=0; jk<32; jk++) counts[jk] = reco.GetZDCScaler(jk);
af6f24c9 1349 fESDZDC->SetZDCScaler(counts);
9e05925b 1350 }
82dffa48 1351
0edb8666 1352 Int_t tdcValues[32][4] = {{0,}};
1353 Float_t tdcCorrected[32][4] = {{0.,}};
32e2fda5 1354 for(Int_t jk=0; jk<32; jk++){
9e05925b 1355 for(Int_t lk=0; lk<4; lk++){
f53e5ecb 1356 tdcValues[jk][lk] = reco.GetZDCTDCData(jk, lk);
27d1ff1f 1357 //Ch debug
1358 //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]);
9e05925b 1359 }
32e2fda5 1360 }
27d1ff1f 1361
1362 // Writing TDC data into ZDC ESDs
44c5ae37 1363 // 4/2/2011 -> Subtracting L0 (tdcValues[15]) instead of ADC gate
1364 // we try to keep the TDC oscillations as low as possible!
f53e5ecb 1365 for(Int_t jk=0; jk<32; jk++){
1366 for(Int_t lk=0; lk<4; lk++){
45b889f8 1367 if(tdcValues[jk][lk]!=0.){
27d1ff1f 1368 tdcCorrected[jk][lk] = 0.025*(tdcValues[jk][lk]-tdcValues[15][0])+fMeanPhase;
1369 // Sep 2011: TDC ch. from 8 to 13 centered around 0 using OCDB
1370 if(jk>=8 && jk<=13) tdcCorrected[jk][lk] = tdcCorrected[jk][lk] - tdcOffset[jk-8];
1371 //Ch. debug
b64f2ac3 1372 //if(jk>=8 && jk<=13) printf(" *** tdcOffset%d %f tdcCorr%d %f \n",jk,tdcOffset[jk-8],tdcCorrected[jk][lk]);
27d1ff1f 1373
1374 }
f53e5ecb 1375 }
1376 }
27d1ff1f 1377
f53e5ecb 1378 fESDZDC->SetZDCTDCData(tdcValues);
1379 fESDZDC->SetZDCTDCCorrected(tdcCorrected);
077d2505 1380 fESDZDC->AliESDZDC::SetBit(AliESDZDC::kCorrectedTDCFilled, reco.GetEnergyFlag());
950abd38 1381 fESDZDC->AliESDZDC::SetBit(AliESDZDC::kEnergyCalibratedSignal, kTRUE);
9e05925b 1382
af6f24c9 1383 if(esd) esd->SetZDCData(fESDZDC);
8309c1ab 1384}
48642b09 1385
1386//_____________________________________________________________________________
78d18275 1387AliCDBStorage* AliZDCReconstructor::SetStorage(const char *uri)
48642b09 1388{
cc2abffd 1389 // Setting the storage
48642b09 1390
78d18275 1391 Bool_t deleteManager = kFALSE;
48642b09 1392
78d18275 1393 AliCDBManager *manager = AliCDBManager::Instance();
1394 AliCDBStorage *defstorage = manager->GetDefaultStorage();
48642b09 1395
78d18275 1396 if(!defstorage || !(defstorage->Contains("ZDC"))){
1397 AliWarning("No default storage set or default storage doesn't contain ZDC!");
1398 manager->SetDefaultStorage(uri);
1399 deleteManager = kTRUE;
1400 }
1401
1402 AliCDBStorage *storage = manager->GetDefaultStorage();
1403
1404 if(deleteManager){
1405 AliCDBManager::Instance()->UnsetDefaultStorage();
1406 defstorage = 0; // the storage is killed by AliCDBManager::Instance()->Destroy()
1407 }
1408
1409 return storage;
1410}
48642b09 1411
78d18275 1412//_____________________________________________________________________________
1e319f71 1413AliZDCPedestals* AliZDCReconstructor::GetPedestalData() const
6024ec85 1414{
1415
1416 // Getting pedestal calibration object for ZDC set
1417
1418 AliCDBEntry *entry = AliCDBManager::Instance()->Get("ZDC/Calib/Pedestals");
f00081c8 1419 if(!entry) AliFatal("No calibration data loaded!");
1420 entry->SetOwner(kFALSE);
6024ec85 1421
1422 AliZDCPedestals *calibdata = dynamic_cast<AliZDCPedestals*> (entry->GetObject());
1423 if(!calibdata) AliFatal("Wrong calibration object in calibration file!");
1424
1425 return calibdata;
1426}
1427
1428//_____________________________________________________________________________
1e319f71 1429AliZDCEnCalib* AliZDCReconstructor::GetEnergyCalibData() const
6024ec85 1430{
1431
1432 // Getting energy and equalization calibration object for ZDC set
1433
dd98e862 1434 AliCDBEntry *entry = AliCDBManager::Instance()->Get("ZDC/Calib/EnergyCalib");
6024ec85 1435 if(!entry) AliFatal("No calibration data loaded!");
f00081c8 1436 entry->SetOwner(kFALSE);
6024ec85 1437
73bc3a3f 1438 AliZDCEnCalib *calibdata = dynamic_cast<AliZDCEnCalib*> (entry->GetObject());
6024ec85 1439 if(!calibdata) AliFatal("Wrong calibration object in calibration file!");
1440
1441 return calibdata;
1442}
1443
73bc3a3f 1444//_____________________________________________________________________________
1e319f71 1445AliZDCTowerCalib* AliZDCReconstructor::GetTowerCalibData() const
73bc3a3f 1446{
1447
1448 // Getting energy and equalization calibration object for ZDC set
1449
dd98e862 1450 AliCDBEntry *entry = AliCDBManager::Instance()->Get("ZDC/Calib/TowerCalib");
73bc3a3f 1451 if(!entry) AliFatal("No calibration data loaded!");
f00081c8 1452 entry->SetOwner(kFALSE);
73bc3a3f 1453
1454 AliZDCTowerCalib *calibdata = dynamic_cast<AliZDCTowerCalib*> (entry->GetObject());
1455 if(!calibdata) AliFatal("Wrong calibration object in calibration file!");
1456
1457 return calibdata;
1458}
1459
1460//_____________________________________________________________________________
0d579f58 1461AliZDCMBCalib* AliZDCReconstructor::GetMBCalibData() const
1462{
1463
1464 // Getting energy and equalization calibration object for ZDC set
1465
1466 AliCDBEntry *entry = AliCDBManager::Instance()->Get("ZDC/Calib/MBCalib");
1467 if(!entry) AliFatal("No calibration data loaded!");
1468 entry->SetOwner(kFALSE);
1469
1470 AliZDCMBCalib *calibdata = dynamic_cast<AliZDCMBCalib*> (entry->GetObject());
1471 if(!calibdata) AliFatal("Wrong calibration object in calibration file!");
1472
1473 return calibdata;
1474}
27d1ff1f 1475
1476//_____________________________________________________________________________
1477AliZDCTDCCalib* AliZDCReconstructor::GetTDCCalibData() const
1478{
1479
1480 // Getting TDC object for ZDC
1481
1482 AliCDBEntry *entry = AliCDBManager::Instance()->Get("ZDC/Calib/TDCCalib");
1483 if(!entry) AliFatal("No calibration data loaded!");
1484 entry->SetOwner(kFALSE);
1485
1486 AliZDCTDCCalib *calibdata = dynamic_cast<AliZDCTDCCalib*> (entry->GetObject());
1487 if(!calibdata) AliFatal("Wrong calibration object in calibration file!");
1488
1489 return calibdata;
1490}