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