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