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