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