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