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