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