1 /**************************************************************************
2 * Copyright(c) 1998-1999, ALICE Experiment at CERN, All rights reserved. *
4 * Author: The ALICE Off-line Project. *
5 * Contributors are mentioned in the code where appropriate. *
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 *
13 * provided "as is" without express or implied warranty. *
14 **************************************************************************/
18 ///////////////////////////////////////////////////////////////////////////////
20 // class for ZDC reconstruction //
22 ///////////////////////////////////////////////////////////////////////////////
27 #include "AliRunLoader.h"
28 #include "AliRawReader.h"
29 #include "AliESDEvent.h"
30 #include "AliESDZDC.h"
31 #include "AliZDCDigit.h"
32 #include "AliZDCRawStream.h"
33 #include "AliZDCReco.h"
34 #include "AliZDCReconstructor.h"
35 #include "AliZDCPedestals.h"
36 #include "AliZDCCalib.h"
37 #include "AliZDCRecParam.h"
40 ClassImp(AliZDCReconstructor)
43 //_____________________________________________________________________________
44 AliZDCReconstructor:: AliZDCReconstructor() :
46 fZNCen(new TF1("fZNCen",
47 "(-2.287920+sqrt(2.287920*2.287920-4*(-0.007629)*(11.921710-x)))/(2*(-0.007629))",0.,164.)),
48 fZNPer(new TF1("fZNPer",
49 "(-37.812280-sqrt(37.812280*37.812280-4*(-0.190932)*(-1709.249672-x)))/(2*(-0.190932))",0.,164.)),
50 fZPCen(new TF1("fZPCen",
51 "(-1.321353+sqrt(1.321353*1.321353-4*(-0.007283)*(3.550697-x)))/(2*(-0.007283))",0.,60.)),
52 fZPPer(new TF1("fZPPer",
53 "(-42.643308-sqrt(42.643308*42.643308-4*(-0.310786)*(-1402.945615-x)))/(2*(-0.310786))",0.,60.)),
54 fZDCCen(new TF1("fZDCCen",
55 "(-1.934991+sqrt(1.934991*1.934991-4*(-0.004080)*(15.111124-x)))/(2*(-0.004080))",0.,225.)),
56 fZDCPer(new TF1("fZDCPer",
57 "(-34.380639-sqrt(34.380639*34.380639-4*(-0.104251)*(-2612.189017-x)))/(2*(-0.104251))",0.,225.)),
58 fbCen(new TF1("fbCen","-0.056923+0.079703*x-0.0004301*x*x+0.000001366*x*x*x",0.,220.)),
59 fbPer(new TF1("fbPer","17.943998-0.046846*x+0.000074*x*x",0.,220.)),
61 fZEMn(new TF1("fZEMn","121.7-0.1934*x+0.00007565*x*x",0.,1200.)),
62 fZEMp(new TF1("fZEMp","80.05-0.1315*x+0.00005327*x*x",0.,1200.)),
63 fZEMsp(new TF1("fZEMsp","201.7-0.325*x+0.0001292*x*x",0.,1200.)),
64 fZEMb(new TF1("fZEMb",
65 "13.83-0.02851*x+5.101e-5*x*x-7.305e-8*x*x*x+5.101e-11*x*x*x*x-1.25e-14*x*x*x*x*x",0.,1200.)),
67 fPedData(GetPedData()),
68 fECalibData(GetECalibData()),
69 fRecParam(GetRecParams())
71 // **** Default constructor
76 //_____________________________________________________________________________
77 AliZDCReconstructor::~AliZDCReconstructor()
97 //_____________________________________________________________________________
98 void AliZDCReconstructor::Reconstruct(TTree* digitsTree, TTree* clustersTree) const
100 // *** Local ZDC reconstruction for digits
101 // Works on the current event
103 // Retrieving calibration data
105 for(Int_t jj=0; jj<48; jj++) meanPed[jj] = fPedData->GetMeanPed(jj);
109 AliZDCDigit* pdigit = &digit;
110 digitsTree->SetBranchAddress("ZDC", &pdigit);
111 //printf("\n\t # of digits in tree: %d\n",(Int_t) digitsTree->GetEntries());
114 Float_t tZN1Corr[10], tZP1Corr[10], tZN2Corr[10], tZP2Corr[10];
115 Float_t dZEM1Corr[2], dZEM2Corr[2], PMRef1[2], PMRef2[2];
116 for(Int_t i=0; i<10; i++){
117 tZN1Corr[i] = tZP1Corr[i] = tZN2Corr[i] = tZP2Corr[i] = 0.;
118 if(i<2) dZEM1Corr[i] = dZEM2Corr[i] = PMRef1[i] = PMRef2[i] = 0.;
121 for (Int_t iDigit = 0; iDigit < (digitsTree->GetEntries()/2); iDigit++) {
122 digitsTree->GetEntry(iDigit);
123 if (!pdigit) continue;
125 Int_t det = digit.GetSector(0);
126 Int_t quad = digit.GetSector(1);
127 Int_t pedindex = -1, kNch = 24;
128 //printf("\n\t Digit #%d det %d quad %d", iDigit, det, quad);
130 if(quad != 5){ // ZDC (not reference PTMs!)
131 if(det == 1){ // *** ZNC
133 tZN1Corr[quad] = (Float_t) (digit.GetADCValue(0)-meanPed[pedindex]);
134 if(tZN1Corr[quad]<0.) tZN1Corr[quad] = 0.;
135 tZN1Corr[quad+5] = (Float_t) (digit.GetADCValue(1)-meanPed[pedindex+kNch]);
136 if(tZN1Corr[quad+5]<0.) tZN1Corr[quad] = 0.;
137 //printf("\t pedindex %d tZN1Corr[%d] = %1.0f tZN1Corr[%d] = %1.0f",
138 // pedindex, quad, tZN1Corr[quad], quad+5, tZN1Corr[quad+5]);
140 else if(det == 2){ // *** ZP1
142 tZP1Corr[quad] = (Float_t) (digit.GetADCValue(0)-meanPed[pedindex]);
143 if(tZP1Corr[quad]<0.) tZP1Corr[quad] = 0.;
144 tZP1Corr[quad+5] = (Float_t) (digit.GetADCValue(1)-meanPed[pedindex+kNch]);
145 if(tZP1Corr[quad+5]<0.) tZP1Corr[quad] = 0.;
146 //printf("\t pedindex %d tZP1Corr[%d] = %1.0f tZP1Corr[%d] = %1.0f",
147 // pedindex, quad, tZP1Corr[quad], quad+5, tZP1Corr[quad+5]);
151 if(quad == 1){ // *** ZEM1
152 dZEM1Corr[0] += (Float_t) (digit.GetADCValue(0)-meanPed[pedindex]);
153 if(dZEM1Corr[0]<0.) dZEM1Corr[0] = 0.;
154 dZEM1Corr[1] += (Float_t) (digit.GetADCValue(1)-meanPed[pedindex+kNch]);
155 if(dZEM1Corr[1]<0.) dZEM1Corr[1] = 0.;
156 //printf("\t pedindex %d tZEM1Corr[%d] = %1.0f tZEM1Corr[%d] = %1.0f",
157 // pedindex, quad, tZEM1Corr[quad], quad+1, tZEM1Corr[quad+1]);
159 else if(quad == 2){ // *** ZEM2
160 dZEM2Corr[0] += (Float_t) (digit.GetADCValue(0)-meanPed[pedindex]);
161 if(dZEM2Corr[0]<0.) dZEM2Corr[0] = 0.;
162 dZEM2Corr[1] += (Float_t) (digit.GetADCValue(1)-meanPed[pedindex+kNch]);
163 if(dZEM2Corr[1]<0.) dZEM2Corr[1] = 0.;
164 //printf("\t pedindex %d tZEM2Corr[%d] = %1.0f tZEM2Corr[%d] = %1.0f",
165 // pedindex, quad, tZEM2Corr[quad], quad+1, tZEM2Corr[quad+1]);
168 else if(det == 4){ // *** ZN2
170 tZN2Corr[quad] = (Float_t) (digit.GetADCValue(0)-meanPed[pedindex]);
171 if(tZN2Corr[quad]<0.) tZN2Corr[quad] = 0.;
172 tZN2Corr[quad+5] = (Float_t) (digit.GetADCValue(1)-meanPed[pedindex+kNch]);
173 if(tZN2Corr[quad+5]<0.) tZN2Corr[quad+5] = 0.;
174 //printf("\t pedindex %d tZN2Corr[%d] = %1.0f tZN2Corr[%d] = %1.0f\n",
175 // pedindex, quad, tZN2Corr[quad], quad+5, tZN2Corr[quad+5]);
177 else if(det == 5){ // *** ZP2
179 tZP2Corr[quad] = (Float_t) (digit.GetADCValue(0)-meanPed[pedindex]);
180 if(tZP2Corr[quad]<0.) tZP2Corr[quad] = 0.;
181 tZP2Corr[quad+5] = (Float_t) (digit.GetADCValue(1)-meanPed[pedindex+kNch]);
182 if(tZP2Corr[quad+5]<0.) tZP2Corr[quad+5] = 0.;
183 //printf("\t pedindex %d tZP2Corr[%d] = %1.0f tZP2Corr[%d] = %1.0f\n",
184 // pedindex, quad, tZP2Corr[quad], quad+5, tZP2Corr[quad+5]);
187 else{ // Reference PMs
188 pedindex = (det-1)/3+22;
190 PMRef1[0] = (Float_t) (digit.GetADCValue(0)-meanPed[pedindex]);
191 if(PMRef1[0]<0.) PMRef1[0] = 0.;
192 PMRef1[1] = (Float_t) (digit.GetADCValue(1)-meanPed[pedindex+kNch]);
193 if(PMRef2[1]<0.) PMRef1[1] = 0.;
196 PMRef2[0] = (Float_t) (digit.GetADCValue(0)-meanPed[pedindex]);
197 if(PMRef2[0]<0.) PMRef2[0] = 0.;
198 PMRef2[1] = (Float_t) (digit.GetADCValue(1)-meanPed[pedindex+kNch]);
199 if(PMRef2[1]<0.) PMRef2[1] = 0.;
204 // reconstruct the event
205 ReconstructEvent(clustersTree, tZN1Corr, tZP1Corr, tZN2Corr, tZP2Corr,
206 dZEM1Corr, dZEM2Corr, PMRef1, PMRef2);
210 //_____________________________________________________________________________
211 void AliZDCReconstructor::Reconstruct(AliRawReader* rawReader, TTree* clustersTree) const
213 // *** ZDC raw data reconstruction
214 // Works on the current event
216 // Retrieving calibration data
218 for(Int_t jj=0; jj<48; jj++) meanPed[jj] = fPedData->GetMeanPed(jj);
222 // loop over raw data
223 Float_t tZN1Corr[10], tZP1Corr[10], tZN2Corr[10], tZP2Corr[10];
224 Float_t dZEM1Corr[2], dZEM2Corr[2], PMRef1[2], PMRef2[2];
225 for(Int_t i=0; i<10; i++){
226 tZN1Corr[i] = tZP1Corr[i] = tZN2Corr[i] = tZP2Corr[i] = 0.;
227 if(i<2) dZEM1Corr[i] = dZEM2Corr[i] = PMRef1[i] = PMRef2[i] = 0.;
230 AliZDCRawStream rawData(rawReader);
232 while (rawData.Next()) {
233 if(rawData.IsADCDataWord()){
234 Int_t det = rawData.GetSector(0);
235 Int_t quad = rawData.GetSector(1);
236 Int_t gain = rawData.GetADCGain();
239 if(quad !=5){ // ZDCs (not reference PTMs)
242 if(gain == 0) tZN1Corr[quad] += (Float_t) (rawData.GetADCValue()-meanPed[pedindex]);
243 else tZN1Corr[quad+5] += (Float_t) (rawData.GetADCValue()-meanPed[pedindex+kNch]);
247 if(gain == 0) tZP1Corr[quad] += (Float_t) (rawData.GetADCValue()-meanPed[pedindex]);
248 else tZP1Corr[quad+5] += (Float_t) (rawData.GetADCValue()-meanPed[pedindex+kNch]);
253 if(gain == 0) dZEM1Corr[0] += (Float_t) (rawData.GetADCValue()-meanPed[pedindex]);
254 else dZEM1Corr[1] += (Float_t) (rawData.GetADCValue()-meanPed[pedindex+kNch]);
257 if(gain == 0) dZEM2Corr[0] += (Float_t) (rawData.GetADCValue()-meanPed[pedindex]);
258 else dZEM2Corr[1] += (Float_t) (rawData.GetADCValue()-meanPed[pedindex+kNch]);
263 if(gain == 0) tZN2Corr[quad] += (Float_t) (rawData.GetADCValue()-meanPed[pedindex]);
264 else tZN2Corr[quad+5] += (Float_t) (rawData.GetADCValue()-meanPed[pedindex+kNch]);
268 if(gain == 0) tZP2Corr[quad] += (Float_t) (rawData.GetADCValue()-meanPed[pedindex]);
269 else tZP2Corr[quad+5] += (Float_t) (rawData.GetADCValue()-meanPed[pedindex+kNch]);
271 //printf("\t AliZDCReconstructor - det %d quad %d res %d -> Ped[%d] = %1.0f\n",
272 // det,quad,gain, pedindex, meanPed[pedindex]);
274 else{ // reference PM
275 pedindex = (det-1)/3 + 22;
277 if(gain==0) PMRef1[0] += (Float_t) (rawData.GetADCValue()-meanPed[pedindex]);
278 else PMRef1[1] += (Float_t) (rawData.GetADCValue()-meanPed[pedindex]);
281 if(gain==0) PMRef2[0] += (Float_t) (rawData.GetADCValue()-meanPed[pedindex]);
282 else PMRef2[1] += (Float_t) (rawData.GetADCValue()-meanPed[pedindex]);
288 // reconstruct the event
289 ReconstructEvent(clustersTree, tZN1Corr, tZP1Corr, tZN2Corr, tZP2Corr,
290 dZEM1Corr, dZEM2Corr, PMRef1, PMRef2);
294 //_____________________________________________________________________________
295 void AliZDCReconstructor::ReconstructEvent(TTree *clustersTree, Float_t* ZN1ADCCorr,
296 Float_t* ZP1ADCCorr, Float_t* ZN2ADCCorr, Float_t* ZP2ADCCorr,
297 Float_t* ZEM1ADCCorr, Float_t* ZEM2ADCCorr, Float_t* PMRef1, Float_t* PMRef2) const
299 // ***** Reconstruct one event
301 // *** RECONSTRUCTION FROM "REAL" DATA
303 // Retrieving calibration data
304 // --- Equalization coefficients ---------------------------------------------
305 Float_t equalCoeffZN1[5], equalCoeffZP1[5], equalCoeffZN2[5], equalCoeffZP2[5];
306 for(Int_t ji=0; ji<5; ji++){
307 equalCoeffZN1[ji] = fECalibData->GetZN1EqualCoeff(ji);
308 equalCoeffZP1[ji] = fECalibData->GetZP1EqualCoeff(ji);
309 equalCoeffZN2[ji] = fECalibData->GetZN2EqualCoeff(ji);
310 equalCoeffZP2[ji] = fECalibData->GetZP2EqualCoeff(ji);
312 // --- Energy calibration factors ------------------------------------
314 for(Int_t ij=0; ij<4; ij++) calibEne[ij] = fECalibData->GetEnCalib(ij);
316 // --- Reconstruction parameters ------------------
317 Float_t endPointZEM = fRecParam->GetZEMEndValue();
318 Float_t cutFractionZEM = fRecParam->GetZEMCutFraction();
319 Float_t dZEMSup = fRecParam->GetDZEMSup();
320 Float_t dZEMInf = fRecParam->GetDZEMInf();
322 Float_t cutValueZEM = endPointZEM*cutFractionZEM;
323 Float_t supValueZEM = cutValueZEM+(endPointZEM*dZEMSup);
324 Float_t infValueZEM = cutValueZEM-(endPointZEM*dZEMInf);
326 Float_t maxValEZN1 = fRecParam->GetEZN1MaxValue();
327 Float_t maxValEZP1 = fRecParam->GetEZP1MaxValue();
328 Float_t maxValEZDC1 = fRecParam->GetEZDC1MaxValue();
329 Float_t maxValEZN2 = fRecParam->GetEZN2MaxValue();
330 Float_t maxValEZP2 = fRecParam->GetEZP2MaxValue();
331 Float_t maxValEZDC2 = fRecParam->GetEZDC2MaxValue();
333 //printf("\n\t AliZDCReconstructor -> ZEMEndPoint %1.0f, ZEMCutValue %1.0f,"
334 // " ZEMSupValue %1.0f, ZEMInfValue %1.0f\n",endPointZEM,cutValueZEM,supValueZEM,infValueZEM);
336 // Equalization of detector responses
337 Float_t equalTowZN1[10], equalTowZN2[10], equalTowZP1[10], equalTowZP2[10];
338 for(Int_t gi=0; gi<5; gi++){
339 equalTowZN1[gi] = ZN1ADCCorr[gi]*equalCoeffZN1[gi];
340 equalTowZN1[gi+5] = ZN1ADCCorr[gi+5]*equalCoeffZN1[gi];
341 equalTowZP1[gi] = ZP1ADCCorr[gi]*equalCoeffZP1[gi];
342 equalTowZP1[gi+5] = ZP1ADCCorr[gi+5]*equalCoeffZP1[gi];
343 equalTowZN2[gi] = ZN2ADCCorr[gi]*equalCoeffZN2[gi];
344 equalTowZN2[gi+5] = ZN2ADCCorr[gi+5]*equalCoeffZN2[gi];
345 equalTowZP2[gi] = ZP2ADCCorr[gi]*equalCoeffZP2[gi];
346 equalTowZP2[gi+5] = ZP2ADCCorr[gi+5]*equalCoeffZP2[gi];
349 // Energy calibration of detector responses
350 Float_t calibTowZN1[10], calibTowZN2[10], calibTowZP1[10], calibTowZP2[10];
351 Float_t calibSumZN1[]={0,0}, calibSumZN2[]={0,0}, calibSumZP1[]={0,0}, calibSumZP2[]={0,0};
352 for(Int_t gi=0; gi<10; gi++){
353 calibTowZN1[gi] = equalTowZN1[gi]*calibEne[0];
354 calibTowZP1[gi] = equalTowZP1[gi]*calibEne[1];
355 calibTowZN2[gi] = equalTowZN2[gi]*calibEne[2];
356 calibTowZP2[gi] = equalTowZP2[gi]*calibEne[3];
359 calibSumZN1[0] += calibTowZN1[gi];
360 calibSumZP1[0] += calibTowZP1[gi];
361 calibSumZN2[0] += calibTowZN2[gi];
362 calibSumZP2[0] += calibTowZP2[gi];
366 calibSumZN1[1] += calibTowZN1[gi];
367 calibSumZP1[1] += calibTowZP1[gi];
368 calibSumZN2[1] += calibTowZN2[gi];
369 calibSumZP2[1] += calibTowZP2[gi];
373 // --- Number of detected spectator nucleons
374 // *** N.B. -> It works only in Pb-Pb
375 Int_t nDetSpecNLeft, nDetSpecPLeft, nDetSpecNRight, nDetSpecPRight;
376 nDetSpecNLeft = (Int_t) (calibSumZN1[0]/2.760);
377 nDetSpecPLeft = (Int_t) (calibSumZP1[0]/2.760);
378 nDetSpecNRight = (Int_t) (calibSumZN2[0]/2.760);
379 nDetSpecPRight = (Int_t) (calibSumZP2[0]/2.760);
380 /*printf("\n\t AliZDCReconstructor -> nDetSpecNLeft %d, nDetSpecPLeft %d,"
381 " nDetSpecNRight %d, nDetSpecPRight %d\n",nDetSpecNLeft, nDetSpecPLeft,
382 nDetSpecNRight, nDetSpecPRight);*/
384 // --- Number of generated spectator nucleons (from HIJING parameterization)
385 Int_t nGenSpecNLeft=0, nGenSpecPLeft=0, nGenSpecLeft=0;
386 Int_t nGenSpecNRight=0, nGenSpecPRight=0, nGenSpecRight=0;
390 Float_t corrADCZEMHG = ZEM1ADCCorr[0] + ZEM2ADCCorr[0];
392 if(corrADCZEMHG > supValueZEM){
393 nGenSpecNLeft = (Int_t) (fZNCen->Eval(calibSumZN1[0]));
394 nGenSpecPLeft = (Int_t) (fZPCen->Eval(calibSumZP1[0]));
395 nGenSpecLeft = (Int_t) (fZDCCen->Eval(calibSumZN1[0]+calibSumZP1[0]));
396 nGenSpecNRight = (Int_t) (fZNCen->Eval(calibSumZN2[0]));
397 nGenSpecPRight = (Int_t) (fZNCen->Eval(calibSumZP2[0]));
398 nGenSpecRight = (Int_t) (fZNCen->Eval(calibSumZN2[0]+calibSumZP2[0]));
399 impPar = fbCen->Eval(calibSumZN1[0]+calibSumZP1[0]);
401 else if(corrADCZEMHG < infValueZEM){
402 nGenSpecNLeft = (Int_t) (fZNPer->Eval(calibSumZN1[0]));
403 nGenSpecPLeft = (Int_t) (fZPPer->Eval(calibSumZP1[0]));
404 nGenSpecLeft = (Int_t) (fZDCPer->Eval(calibSumZN1[0]+calibSumZP1[0]));
405 impPar = fbPer->Eval(calibSumZN1[0]+calibSumZP1[0]);
407 else if(corrADCZEMHG >= infValueZEM && corrADCZEMHG <= supValueZEM){
408 nGenSpecNLeft = (Int_t) (fZEMn->Eval(corrADCZEMHG));
409 nGenSpecPLeft = (Int_t) (fZEMp->Eval(corrADCZEMHG));
410 nGenSpecLeft = (Int_t)(fZEMsp->Eval(corrADCZEMHG));
411 impPar = fZEMb->Eval(corrADCZEMHG);
414 if(calibSumZN1[0]/maxValEZN1>1.) nGenSpecNLeft = (Int_t) (fZEMn->Eval(corrADCZEMHG));
415 if(calibSumZP1[0]/maxValEZP1>1.) nGenSpecPLeft = (Int_t) (fZEMp->Eval(corrADCZEMHG));
416 if((calibSumZN1[0]+calibSumZP1[0]/maxValEZDC1)>1.){
417 nGenSpecLeft = (Int_t)(fZEMsp->Eval(corrADCZEMHG));
418 impPar = fZEMb->Eval(corrADCZEMHG);
420 if(calibSumZN2[0]/maxValEZN2>1.) nGenSpecNRight = (Int_t) (fZEMn->Eval(corrADCZEMHG));
421 if(calibSumZP2[0]/maxValEZP2>1.) nGenSpecPRight = (Int_t) (fZEMp->Eval(corrADCZEMHG));
422 if((calibSumZN2[0]+calibSumZP2[0]/maxValEZDC2)>1.) nGenSpecRight = (Int_t)(fZEMsp->Eval(corrADCZEMHG));
424 if(nGenSpecNLeft>125) nGenSpecNLeft=125;
425 else if(nGenSpecNLeft<0) nGenSpecNLeft=0;
426 if(nGenSpecPLeft>82) nGenSpecPLeft=82;
427 else if(nGenSpecPLeft<0) nGenSpecPLeft=0;
428 if(nGenSpecLeft>207) nGenSpecLeft=207;
429 else if(nGenSpecLeft<0) nGenSpecLeft=0;
431 // --- Number of generated participants (from HIJING parameterization)
432 Int_t nPart, nPartTotLeft, nPartTotRight;
433 nPart = 207-nGenSpecNLeft-nGenSpecPLeft;
434 nPartTotLeft = 207-nGenSpecLeft;
435 nPartTotRight = 207-nGenSpecRight;
437 if(nPartTotLeft<0) nPartTotLeft=0;
438 if(nPartTotRight<0) nPartTotRight=0;
441 /*printf("\n\t AliZDCReconstructor -> calibSumZN1[0] %1.0f, calibSumZP1[0] %1.0f,"
442 " calibSumZN2[0] %1.0f, calibSumZP2[0] %1.0f, corrADCZEMHG %1.0f\n",
443 calibSumZN1[0],calibSumZP1[0],calibSumZN2[0],calibSumZP2[0],corrADCZEMHG);
444 printf("\t AliZDCReconstructor -> nGenSpecNLeft %d, nGenSpecPLeft %d, nGenSpecLeft %d\n"
445 "\t\t nGenSpecNRight %d, nGenSpecPRight %d, nGenSpecRight %d\n",
446 nGenSpecNLeft, nGenSpecPLeft, nGenSpecLeft,
447 nGenSpecNRight, nGenSpecPRight, nGenSpecRight);
448 printf("\t AliZDCReconstructor -> NpartL %d, NpartR %d, b %1.2f fm\n\n",nPartTotLeft, nPartTotRight, impPar);
451 // create the output tree
452 AliZDCReco reco(calibSumZN1, calibSumZP1, calibSumZN2, calibSumZP2,
453 calibTowZN1, calibTowZN2, calibTowZP1, calibTowZP2,
454 ZEM1ADCCorr, ZEM2ADCCorr, PMRef1, PMRef2,
455 nDetSpecNLeft, nDetSpecPLeft, nDetSpecNRight, nDetSpecPRight,
456 nGenSpecNLeft, nGenSpecPLeft, nGenSpecLeft, nGenSpecNRight,
457 nGenSpecPRight, nGenSpecRight, nPartTotLeft, nPartTotRight, impPar);
459 AliZDCReco* preco = &reco;
460 const Int_t kBufferSize = 4000;
461 clustersTree->Branch("ZDC", "AliZDCReco", &preco, kBufferSize);
463 // write the output tree
464 clustersTree->Fill();
467 //_____________________________________________________________________________
468 void AliZDCReconstructor::FillZDCintoESD(TTree *clustersTree, AliESDEvent* esd) const
470 // fill energies and number of participants to the ESD
473 AliZDCReco* preco = &reco;
474 clustersTree->SetBranchAddress("ZDC", &preco);
476 clustersTree->GetEntry(0);
478 AliESDZDC * esdzdc = esd->GetESDZDC();
479 Float_t tZN1Ene[5], tZN2Ene[5], tZP1Ene[5], tZP2Ene[5];
480 Float_t tZN1EneLR[5], tZN2EneLR[5], tZP1EneLR[5], tZP2EneLR[5];
481 for(Int_t i=0; i<5; i++){
482 tZN1Ene[i] = reco.GetZN1HREnTow(i);
483 tZN2Ene[i] = reco.GetZN2HREnTow(i);
484 tZP1Ene[i] = reco.GetZP1HREnTow(i);
485 tZP2Ene[i] = reco.GetZP2HREnTow(i);
487 tZN1EneLR[i] = reco.GetZN1LREnTow(i);
488 tZN2EneLR[i] = reco.GetZN2LREnTow(i);
489 tZP1EneLR[i] = reco.GetZP1LREnTow(i);
490 tZP2EneLR[i] = reco.GetZP2LREnTow(i);
492 esdzdc->SetZN1TowerEnergy(tZN1Ene);
493 esdzdc->SetZN2TowerEnergy(tZN2Ene);
494 esdzdc->SetZP1TowerEnergy(tZP1Ene);
495 esdzdc->SetZP2TowerEnergy(tZP2Ene);
496 esdzdc->SetZN1TowerEnergyLR(tZN1EneLR);
497 esdzdc->SetZN2TowerEnergyLR(tZN2EneLR);
498 esdzdc->SetZP1TowerEnergyLR(tZP1EneLR);
499 esdzdc->SetZP2TowerEnergyLR(tZP2EneLR);
501 esd->SetZDC(reco.GetZN1HREnergy(), reco.GetZP1HREnergy(), reco.GetZEM1HRsignal(),
502 reco.GetZEM2HRsignal(), reco.GetZN2HREnergy(), reco.GetZP2HREnergy(),
503 reco.GetNPartLeft());
508 //_____________________________________________________________________________
509 AliCDBStorage* AliZDCReconstructor::SetStorage(const char *uri)
511 // Setting the storage
513 Bool_t deleteManager = kFALSE;
515 AliCDBManager *manager = AliCDBManager::Instance();
516 AliCDBStorage *defstorage = manager->GetDefaultStorage();
518 if(!defstorage || !(defstorage->Contains("ZDC"))){
519 AliWarning("No default storage set or default storage doesn't contain ZDC!");
520 manager->SetDefaultStorage(uri);
521 deleteManager = kTRUE;
524 AliCDBStorage *storage = manager->GetDefaultStorage();
527 AliCDBManager::Instance()->UnsetDefaultStorage();
528 defstorage = 0; // the storage is killed by AliCDBManager::Instance()->Destroy()
534 //_____________________________________________________________________________
535 AliZDCPedestals* AliZDCReconstructor::GetPedData() const
538 // Getting pedestal calibration object for ZDC set
540 AliCDBEntry *entry = AliCDBManager::Instance()->Get("ZDC/Calib/Pedestals");
541 if(!entry) AliFatal("No calibration data loaded!");
543 AliZDCPedestals *calibdata = dynamic_cast<AliZDCPedestals*> (entry->GetObject());
544 if(!calibdata) AliFatal("Wrong calibration object in calibration file!");
549 //_____________________________________________________________________________
550 AliZDCCalib* AliZDCReconstructor::GetECalibData() const
553 // Getting energy and equalization calibration object for ZDC set
555 AliCDBEntry *entry = AliCDBManager::Instance()->Get("ZDC/Calib/Calib");
556 if(!entry) AliFatal("No calibration data loaded!");
558 AliZDCCalib *calibdata = dynamic_cast<AliZDCCalib*> (entry->GetObject());
559 if(!calibdata) AliFatal("Wrong calibration object in calibration file!");
564 //_____________________________________________________________________________
565 AliZDCRecParam* AliZDCReconstructor::GetRecParams() const
568 // Getting energy and equalization calibration object for ZDC set
570 AliCDBEntry *entry = AliCDBManager::Instance()->Get("ZDC/Calib/RecParam");
571 if(!entry) AliFatal("No calibration data loaded!");
573 AliZDCRecParam *calibdata = dynamic_cast<AliZDCRecParam*> (entry->GetObject());
574 if(!calibdata) AliFatal("Wrong calibration object in calibration file!");