small fix (Giuseppe)
[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// //
20// class for ZDC reconstruction //
21// //
22///////////////////////////////////////////////////////////////////////////////
23
24
25#include <TF1.h>
26
27#include "AliRunLoader.h"
28#include "AliRawReader.h"
29#include "AliESDEvent.h"
30#include "AliZDCDigit.h"
31#include "AliZDCRawStream.h"
32#include "AliZDCReco.h"
33#include "AliZDCReconstructor.h"
6024ec85 34#include "AliZDCPedestals.h"
35#include "AliZDCCalib.h"
36#include "AliZDCRecParam.h"
f5d41205 37
38
39ClassImp(AliZDCReconstructor)
40
41
42//_____________________________________________________________________________
43AliZDCReconstructor:: AliZDCReconstructor() :
44
45 fZNCen(new TF1("fZNCen",
46 "(-2.287920+sqrt(2.287920*2.287920-4*(-0.007629)*(11.921710-x)))/(2*(-0.007629))",0.,164.)),
47 fZNPer(new TF1("fZNPer",
48 "(-37.812280-sqrt(37.812280*37.812280-4*(-0.190932)*(-1709.249672-x)))/(2*(-0.190932))",0.,164.)),
49 fZPCen(new TF1("fZPCen",
50 "(-1.321353+sqrt(1.321353*1.321353-4*(-0.007283)*(3.550697-x)))/(2*(-0.007283))",0.,60.)),
51 fZPPer(new TF1("fZPPer",
52 "(-42.643308-sqrt(42.643308*42.643308-4*(-0.310786)*(-1402.945615-x)))/(2*(-0.310786))",0.,60.)),
53 fZDCCen(new TF1("fZDCCen",
54 "(-1.934991+sqrt(1.934991*1.934991-4*(-0.004080)*(15.111124-x)))/(2*(-0.004080))",0.,225.)),
55 fZDCPer(new TF1("fZDCPer",
56 "(-34.380639-sqrt(34.380639*34.380639-4*(-0.104251)*(-2612.189017-x)))/(2*(-0.104251))",0.,225.)),
57 fbCen(new TF1("fbCen","-0.056923+0.079703*x-0.0004301*x*x+0.000001366*x*x*x",0.,220.)),
58 fbPer(new TF1("fbPer","17.943998-0.046846*x+0.000074*x*x",0.,220.)),
59 //
60 fZEMn(new TF1("fZEMn","121.7-0.1934*x+0.00007565*x*x",0.,1200.)),
61 fZEMp(new TF1("fZEMp","80.05-0.1315*x+0.00005327*x*x",0.,1200.)),
62 fZEMsp(new TF1("fZEMsp","201.7-0.325*x+0.0001292*x*x",0.,1200.)),
63 fZEMb(new TF1("fZEMb",
64 "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.)),
65 //
6024ec85 66 fPedData(GetPedData()),
67 fECalibData(GetECalibData()),
68 fRecParam(GetRecParams())
f5d41205 69{
70 // **** Default constructor
71
72}
73
74
75//_____________________________________________________________________________
76AliZDCReconstructor::~AliZDCReconstructor()
77{
78// destructor
79
80 delete fZNCen;
81 delete fZNPer;
82 delete fZPCen;
83 delete fZPPer;
84 delete fZDCCen;
85 delete fZDCPer;
86 delete fbCen;
87 delete fbPer;
88 delete fZEMn;
89 delete fZEMp;
90 delete fZEMsp;
91 delete fZEMb;
92
93}
94
95
96//_____________________________________________________________________________
97void AliZDCReconstructor::Reconstruct(TTree* digitsTree, TTree* clustersTree) const
98{
99 // *** Local ZDC reconstruction for digits
100 // Works on the current event
101
102 // Retrieving calibration data
103 Float_t meanPed[47];
6024ec85 104 for(Int_t jj=0; jj<47; jj++) meanPed[jj] = fPedData->GetMeanPed(jj);
f5d41205 105
106 // get digits
107 AliZDCDigit digit;
108 AliZDCDigit* pdigit = &digit;
109 digitsTree->SetBranchAddress("ZDC", &pdigit);
110
111 // loop over digits
112 Float_t tZN1CorrHG[]={0.,0.,0.,0.,0.}, tZP1CorrHG[]={0.,0.,0.,0.,0.};
113 Float_t dZEMCorrHG=0.;
114 Float_t tZN2CorrHG[]={0.,0.,0.,0.,0.}, tZP2CorrHG[]={0.,0.,0.,0.,0.};
115 Float_t tZN1CorrLG[]={0.,0.,0.,0.,0.}, tZP1CorrLG[]={0.,0.,0.,0.,0.};
116 Float_t dZEMCorrLG=0.;
117 Float_t tZN2CorrLG[]={0.,0.,0.,0.,0.}, tZP2CorrLG[]={0.,0.,0.,0.,0.};
118
119 //printf("\n\t # of digits in tree: %d\n",(Int_t) digitsTree->GetEntries());
120 for (Int_t iDigit = 0; iDigit < (digitsTree->GetEntries()/2); iDigit++) {
121 digitsTree->GetEntry(iDigit);
122 if (!pdigit) continue;
123 //pdigit->Print("");
124 //
125 Int_t det = digit.GetSector(0);
126 Int_t quad = digit.GetSector(1);
127 Int_t pedindex = -1;
128 //printf("\n\t #%d det %d quad %d", iDigit, det, quad);
129 //
130 if(det == 1){ // *** ZN1
131 pedindex = quad;
132 tZN1CorrHG[quad] = (Float_t) (digit.GetADCValue(0)-meanPed[pedindex]);
133 if(tZN1CorrHG[quad]<0.) tZN1CorrHG[quad] = 0.;
134 tZN1CorrLG[quad] = (Float_t) (digit.GetADCValue(1)-meanPed[pedindex+5]);
135 if(tZN1CorrLG[quad]<0.) tZN1CorrLG[quad] = 0.;
136 //printf("\t pedindex %d tZN1CorrHG[%d] = %1.0f tZN1CorrLG[%d] = %1.0f",
137 // pedindex, quad, tZN1CorrHG[quad], quad, tZN1CorrLG[quad]);
138 }
139 else if(det == 2){ // *** ZP1
140 pedindex = quad+10;
141 tZP1CorrHG[quad] = (Float_t) (digit.GetADCValue(0)-meanPed[pedindex]);
142 if(tZP1CorrLG[quad]<0.) tZP1CorrLG[quad] = 0.;
143 tZP1CorrLG[quad] = (Float_t) (digit.GetADCValue(1)-meanPed[pedindex+5]);
144 if(tZP1CorrHG[quad]<0.) tZP1CorrHG[quad] = 0.;
145 //printf("\t pedindex %d tZP1CorrHG[%d] = %1.0f tZP1CorrLG[%d] = %1.0f",
146 // pedindex, quad, tZP1CorrHG[quad], quad, tZP1CorrLG[quad]);
147 }
148 else if(det == 3){
149 if(quad == 1){ // *** ZEM1
150 pedindex = quad+19;
151 dZEMCorrHG += (Float_t) (digit.GetADCValue(0)-meanPed[pedindex]);
152 if(dZEMCorrHG<0.) dZEMCorrHG = 0.;
153 dZEMCorrLG += (Float_t) (digit.GetADCValue(1)-meanPed[pedindex+2]);
154 if(dZEMCorrLG<0.) dZEMCorrLG = 0.;
155 //printf("\t pedindex %d ADC(0) = %d ped = %1.0f ADCCorr = %1.0f\n",
156 // pedindex, digit.GetADCValue(0), meanPed[pedindex], dZEMCorrHG);
157 //printf("\t pedindex %d ADC(1) = %d ped = %1.0f ADCCorr = %1.0f\n",
158 // pedindex+2, digit.GetADCValue(1), meanPed[pedindex+2], dZEMCorrLG);
159 ////printf("\t pedindex %d dZEMCorrHG = %1.0f dZEMCorrLG = %1.0f\n", pedindex, dZEMCorrHG, dZEMCorrLG);
160 }
161 else if(quad == 2){ // *** ZEM1
162 pedindex = quad+19;
163 dZEMCorrHG += (Float_t) (digit.GetADCValue(0)-meanPed[pedindex]);
164 if(dZEMCorrHG<0.) dZEMCorrHG = 0.;
165 dZEMCorrLG += (Float_t) (digit.GetADCValue(1)-meanPed[pedindex+2]);
166 if(dZEMCorrLG<0.) dZEMCorrLG = 0.;
167 //printf("\t pedindex %d ADC(0) = %d ped = %1.0f ADCCorr = %1.0f\n",
168 // pedindex, digit.GetADCValue(0), meanPed[pedindex], dZEMCorrHG);
169 //printf("\t pedindex %d ADC(1) = %d ped = %1.0f ADCCorr = %1.0f\n",
170 // pedindex+2, digit.GetADCValue(1),meanPed[pedindex+2], dZEMCorrLG);
171 ////printf("\t pedindex %d dZEMCorrHG = %1.0f dZEMCorrLG = %1.0f\n", pedindex, dZEMCorrHG, dZEMCorrLG);
172 }
173 }
174 else if(det == 4){ // *** ZN2
175 pedindex = quad+24;
176 tZN2CorrHG[quad] = (Float_t) (digit.GetADCValue(0)-meanPed[pedindex]);
177 if(tZN2CorrHG[quad]<0.) tZN2CorrHG[quad] = 0.;
178 tZN2CorrLG[quad] = (Float_t) (digit.GetADCValue(1)-meanPed[pedindex+5]);
179 if(tZN2CorrLG[quad]<0.) tZN2CorrLG[quad] = 0.;
180 //printf("\t pedindex %d tZN2CorrHG[%d] = %1.0f tZN2CorrLG[%d] = %1.0f\n",
181 // pedindex, quad, tZN2CorrHG[quad], quad, tZN2CorrLG[quad]);
182 }
183 else if(det == 5){ // *** ZP2
184 pedindex = quad+34;
185 tZP2CorrHG[quad] = (Float_t) (digit.GetADCValue(0)-meanPed[pedindex]);
186 if(tZP2CorrHG[quad]<0.) tZP2CorrHG[quad] = 0.;
187 tZP2CorrLG[quad] = (Float_t) (digit.GetADCValue(1)-meanPed[pedindex+5]);
188 if(tZP2CorrLG[quad]<0.) tZP2CorrLG[quad] = 0.;
189 //printf("\t pedindex %d tZP2CorrHG[%d] = %1.0f tZP2CorrLG[%d] = %1.0f\n",
190 // pedindex, quad, tZP2CorrHG[quad], quad, tZP2CorrLG[quad]);
191 }
192 }
193
194 // reconstruct the event
195 ReconstructEvent(clustersTree, tZN1CorrHG, tZP1CorrHG, tZN2CorrHG,
196 tZP2CorrHG, tZN1CorrLG, tZP1CorrLG, tZN2CorrLG,
197 tZP2CorrLG, dZEMCorrHG);
198
199}
200
201//_____________________________________________________________________________
202void AliZDCReconstructor::Reconstruct(AliRawReader* rawReader, TTree* clustersTree) const
203{
204 // *** ZDC raw data reconstruction
205 // Works on the current event
206
207 // Retrieving calibration data
208 Float_t meanPed[47];
6024ec85 209 for(Int_t jj=0; jj<47; jj++) meanPed[jj] = fPedData->GetMeanPed(jj);
f5d41205 210
211 rawReader->Reset();
212
213 // loop over raw data rawDatas
214 Float_t tZN1CorrHG[]={0.,0.,0.,0.,0.}, tZP1CorrHG[]={0.,0.,0.,0.,0.};
215 Float_t dZEMCorrHG=0.;
216 Float_t tZN2CorrHG[]={0.,0.,0.,0.,0.}, tZP2CorrHG[]={0.,0.,0.,0.,0.};
217 Float_t tZN1CorrLG[]={0.,0.,0.,0.,0.}, tZP1CorrLG[]={0.,0.,0.,0.,0.};
218 Float_t dZEMCorrLG=0.;
219 Float_t tZN2CorrLG[]={0.,0.,0.,0.,0.}, tZP2CorrLG[]={0.,0.,0.,0.,0.};
220 //
221 AliZDCRawStream rawData(rawReader);
222 while (rawData.Next()) {
223 if(rawData.IsADCDataWord()){
224 Int_t det = rawData.GetSector(0);
225 Int_t quad = rawData.GetSector(1);
226 Int_t gain = rawData.GetADCGain();
227 Int_t pedindex;
228 //
229 if(det == 1){
230 pedindex = quad;
231 if(gain == 0) tZN1CorrHG[quad] += (Float_t) (rawData.GetADCValue()-meanPed[pedindex]);
232 else tZN1CorrLG[quad] += (Float_t) (rawData.GetADCValue()-meanPed[pedindex+5]);
233 }
234 else if(det == 2){
235 pedindex = quad+10;
236 if(gain == 0) tZP1CorrHG[quad] += (Float_t) (rawData.GetADCValue()-meanPed[pedindex]);
237 else tZP1CorrLG[quad] += (Float_t) (rawData.GetADCValue()-meanPed[pedindex+5]);
238 }
239 else if(det == 3){
240 if(quad==1){
241 pedindex = quad+20;
242 if(gain == 0) dZEMCorrHG += (Float_t) (rawData.GetADCValue()-meanPed[pedindex]);
243 else dZEMCorrLG += (Float_t) (rawData.GetADCValue()-meanPed[pedindex+2]);
244 }
245 else if(quad==2){
246 pedindex = rawData.GetSector(1)+21;
247 if(gain == 0) dZEMCorrHG += (Float_t) (rawData.GetADCValue()-meanPed[pedindex]);
248 else dZEMCorrLG += (Float_t) (rawData.GetADCValue()-meanPed[pedindex+2]);
249 }
250 }
251 else if(det == 4){
252 pedindex = rawData.GetSector(1)+24;
253 if(gain == 0) tZN2CorrHG[quad] += (Float_t) (rawData.GetADCValue()-meanPed[pedindex]);
254 else tZN2CorrLG[quad] += (Float_t) (rawData.GetADCValue()-meanPed[pedindex+2]);
255 }
256 else if(det == 5){
257 pedindex = rawData.GetSector(1)+34;
258 if(gain == 0) tZP2CorrHG[quad] += (Float_t) (rawData.GetADCValue()-meanPed[pedindex]);
259 else tZP2CorrLG[quad] += (Float_t) (rawData.GetADCValue()-meanPed[pedindex+5]);
260 }
261 }
262 }
263
264 // reconstruct the event
265 ReconstructEvent(clustersTree, tZN1CorrHG, tZP1CorrHG, tZN2CorrHG,
266 tZP2CorrHG, tZN1CorrLG, tZP1CorrLG, tZN2CorrLG,
267 tZP2CorrLG, dZEMCorrHG);
268
269}
270
271//_____________________________________________________________________________
272void AliZDCReconstructor::ReconstructEvent(TTree *clustersTree,
273 Float_t* ZN1ADCCorrHG, Float_t* ZP1ADCCorrHG,
274 Float_t* ZN2ADCCorrHG, Float_t* ZP2ADCCorrHG,
275 Float_t* ZN1ADCCorrLG, Float_t* ZP1ADCCorrLG,
276 Float_t* ZN2ADCCorrLG, Float_t* ZP2ADCCorrLG,
277 Float_t corrADCZEMHG) const
278{
279 // ***** Reconstruct one event
280
281 // *** RECONSTRUCTION FROM SIMULATED DATA
282 // It passes trhough the no. of phe which is known from simulations
283 // --- ADCchannel -> photoelectrons
284 // NB-> PM gain = 10^(5), ADC resolution = 6.4*10^(-7)
285 // Move to V965 (E.S.,15/09/04) NB-> PM gain = 10^(5), ADC resolution = 8*10^(-7)
286 //Float_t zn1phe, zp1phe, zemphe, zn2phe, zp2phe, convFactor = 0.08;
287 //zn1phe = ZN1Corr/convFactor;
288 //zp1phe = ZP1Corr/convFactor;
289 //zemphe = ZEMCorr/convFactor;
290 //zn2phe = ZN2Corr/convFactor;
291 //zp2phe = ZP2Corr/convFactor;
292 ////if AliDebug(1,Form("\n znphe = %f, zpphe = %f, zemphe = %f\n",znphe, zpphe, zemphe);
293 //
294 //// --- Energy calibration
295 //// Conversion factors for hadronic ZDCs goes from phe yield to TRUE
296 //// incident energy (conversion from GeV to TeV is included); while for EM
297 //// calos conversion is from light yield to detected energy calculated by
298 //// GEANT NB -> ZN and ZP conversion factors are constant since incident
299 //// spectators have all the same energy, ZEM energy is obtained through a
300 //// fit over the whole range of incident particle energies
301 //// (obtained with full HIJING simulations)
302 //Float_t zn1energy, zp1energy, zemenergy, zdc1energy, zn2energy, zp2energy, zdc2energy;
303 //Float_t zn1phexTeV=329., zp1phexTeV=369., zn2phexTeV=329., zp2phexTeV=369.;
304 //zn1energy = zn1phe/zn1phexTeV;
305 //zp1energy = zp1phe/zp1phexTeV;
306 //zdc1energy = zn1energy+zp1energy;
307 //zn2energy = zn2phe/zn2phexTeV;
308 //zp2energy = zp2phe/zp2phexTeV;
309 //zdc2energy = zn2energy+zp2energy;
310 //zemenergy = -4.81+0.3238*zemphe;
311 //if(zemenergy<0) zemenergy=0;
312 //// if AliDebug(1,Form(" znenergy = %f TeV, zpenergy = %f TeV, zdcenergy = %f GeV, "
313 //// "\n zemenergy = %f TeV\n", znenergy, zpenergy,
314 //// zdcenergy, zemenergy);
315 //// if(zdcenergy==0)
316 //// if AliDebug(1,Form("\n\n ### ATTENZIONE!!! -> ev# %d: znenergy = %f TeV, zpenergy = %f TeV, zdcenergy = %f GeV, "
317 //// " zemenergy = %f TeV\n\n", fMerger->EvNum(), znenergy, zpenergy, zdcenergy, zemenergy);
318
319 //
320 // *** RECONSTRUCTION FROM "REAL" DATA
321 //
322 // Retrieving calibration data
84d6255e 323 // --- Equalization coefficients ---------------------------------------------
f5d41205 324 Float_t equalCoeffZN1[5], equalCoeffZP1[5], equalCoeffZN2[5], equalCoeffZP2[5];
325 for(Int_t ji=0; ji<5; ji++){
6024ec85 326 equalCoeffZN1[ji] = fECalibData->GetZN1EqualCoeff(ji);
327 equalCoeffZP1[ji] = fECalibData->GetZP1EqualCoeff(ji);
328 equalCoeffZN2[ji] = fECalibData->GetZN2EqualCoeff(ji);
329 equalCoeffZP2[ji] = fECalibData->GetZP2EqualCoeff(ji);
f5d41205 330 }
84d6255e 331 // --- Energy calibration factors ------------------------------------
f5d41205 332 Float_t calibEne[4];
6024ec85 333 for(Int_t ij=0; ij<4; ij++) calibEne[ij] = fECalibData->GetEnCalib(ij);
f5d41205 334 //
84d6255e 335 // --- Reconstruction parameters ------------------
6024ec85 336 Float_t endPointZEM = fRecParam->GetZEMEndValue();
337 Float_t cutFractionZEM = fRecParam->GetZEMCutFraction();
338 Float_t dZEMSup = fRecParam->GetDZEMSup();
339 Float_t dZEMInf = fRecParam->GetDZEMInf();
f5d41205 340 //
341 Float_t cutValueZEM = endPointZEM*cutFractionZEM;
342 Float_t supValueZEM = cutValueZEM+(endPointZEM*dZEMSup);
343 Float_t infValueZEM = cutValueZEM-(endPointZEM*dZEMInf);
344 //
6024ec85 345 Float_t maxValEZN1 = fRecParam->GetEZN1MaxValue();
346 Float_t maxValEZP1 = fRecParam->GetEZP1MaxValue();
347 Float_t maxValEZDC1 = fRecParam->GetEZDC1MaxValue();
348 Float_t maxValEZN2 = fRecParam->GetEZN2MaxValue();
349 Float_t maxValEZP2 = fRecParam->GetEZP2MaxValue();
350 Float_t maxValEZDC2 = fRecParam->GetEZDC2MaxValue();
f5d41205 351 //
352 //printf("\n\t AliZDCReconstructor -> ZEMEndPoint %1.0f, ZEMCutValue %1.0f,"
353 // " ZEMSupValue %1.0f, ZEMInfValue %1.0f\n",endPointZEM,cutValueZEM,supValueZEM,infValueZEM);
354
355 // Equalization of detector responses
356 Float_t equalTowZN1HG[5], equalTowZN2HG[5], equalTowZP1HG[5], equalTowZP2HG[5];
357 Float_t equalTowZN1LG[5], equalTowZN2LG[5], equalTowZP1LG[5], equalTowZP2LG[5];
358 for(Int_t gi=0; gi<5; gi++){
359 equalTowZN1HG[gi] = ZN1ADCCorrHG[gi]*equalCoeffZN1[gi];
360 equalTowZP1HG[gi] = ZP1ADCCorrHG[gi]*equalCoeffZP1[gi];
361 equalTowZN2HG[gi] = ZN2ADCCorrHG[gi]*equalCoeffZN2[gi];
362 equalTowZP2HG[gi] = ZP2ADCCorrHG[gi]*equalCoeffZP2[gi];
363 //
364 equalTowZN1LG[gi] = ZN1ADCCorrLG[gi]*equalCoeffZN1[gi];
365 equalTowZP1LG[gi] = ZP1ADCCorrLG[gi]*equalCoeffZP1[gi];
366 equalTowZN2LG[gi] = ZN2ADCCorrLG[gi]*equalCoeffZN2[gi];
367 equalTowZP2LG[gi] = ZP2ADCCorrLG[gi]*equalCoeffZP2[gi];
368 }
369
370 // Energy calibration of detector responses
371 Float_t calibTowZN1HG[5], calibTowZN2HG[5], calibTowZP1HG[5], calibTowZP2HG[5];
372 Float_t calibSumZN1HG=0., calibSumZN2HG=0., calibSumZP1HG=0., calibSumZP2HG=0.;
373 Float_t calibTowZN1LG[5], calibTowZN2LG[5], calibTowZP1LG[5], calibTowZP2LG[5];
374 Float_t calibSumZN1LG=0., calibSumZN2LG=0., calibSumZ12LG=0., calibSumZP2LG=0.;
375 for(Int_t gi=0; gi<5; gi++){
376 calibTowZN1HG[gi] = equalTowZN1HG[gi]*calibEne[0];
377 calibTowZP1HG[gi] = equalTowZP1HG[gi]*calibEne[1];
378 calibTowZN2HG[gi] = equalTowZN2HG[gi]*calibEne[2];
379 calibTowZP2HG[gi] = equalTowZP2HG[gi]*calibEne[3];
380 calibSumZN1HG += calibTowZN1HG[gi];
381 calibSumZP1HG += calibTowZP1HG[gi];
382 calibSumZN2HG += calibTowZN2HG[gi];
383 calibSumZP2HG += calibTowZP2HG[gi];
384 //
385 calibTowZN1LG[gi] = equalTowZN1LG[gi]*calibEne[0];
386 calibTowZP1LG[gi] = equalTowZP1LG[gi]*calibEne[1];
387 calibTowZN2LG[gi] = equalTowZN2LG[gi]*calibEne[2];
388 calibTowZP2LG[gi] = equalTowZP2LG[gi]*calibEne[3];
389 calibSumZN1LG += calibTowZN1LG[gi];
390 calibSumZ12LG += calibTowZP1LG[gi];
391 calibSumZN2LG += calibTowZN2LG[gi];
392 calibSumZP2LG += calibTowZP2LG[gi];
393 }
394
395 // --- Number of detected spectator nucleons
396 // *** N.B. -> It works only in Pb-Pb
397 Int_t nDetSpecNLeft, nDetSpecPLeft, nDetSpecNRight, nDetSpecPRight;
398 nDetSpecNLeft = (Int_t) (calibSumZN1HG/2.760);
399 nDetSpecPLeft = (Int_t) (calibSumZP1HG/2.760);
400 nDetSpecNRight = (Int_t) (calibSumZN2HG/2.760);
401 nDetSpecPRight = (Int_t) (calibSumZP2HG/2.760);
402 /*printf("\n\t AliZDCReconstructor -> nDetSpecNLeft %d, nDetSpecPLeft %d,"
403 " nDetSpecNRight %d, nDetSpecPRight %d\n",nDetSpecNLeft, nDetSpecPLeft,
404 nDetSpecNRight, nDetSpecPRight);*/
405
406 // --- Number of generated spectator nucleons (from HIJING parameterization)
407 Int_t nGenSpecNLeft=0, nGenSpecPLeft=0, nGenSpecLeft=0;
408 Int_t nGenSpecNRight=0, nGenSpecPRight=0, nGenSpecRight=0;
409 Double_t impPar=0.;
410 //
411 // *** RECONSTRUCTION FROM SIMULATED DATA
412 // Cut value for Ezem (GeV)
413 // ### Results from production -> 0<b<18 fm (Apr 2002)
414 /*Float_t eZEMCut = 420.;
415 Float_t deltaEZEMSup = 690.;
416 Float_t deltaEZEMInf = 270.;
417 if(zemenergy > (eZEMCut+deltaEZEMSup)){
418 nGenSpecNLeft = (Int_t) (fZNCen->Eval(ZN1CalibSum));
419 nGenSpecPLeft = (Int_t) (fZPCen->Eval(ZP1CalibSum));
420 nGenSpecLeft = (Int_t) (fZDCCen->Eval(ZN1CalibSum+ZP1CalibSum));
421 nGenSpecNRight = (Int_t) (fZNCen->Eval(ZN2CalibSum));
422 nGenSpecPRight = (Int_t) (fZNCen->Eval(ZP2CalibSum));
423 nGenSpecRight = (Int_t) (fZNCen->Eval(ZN2CalibSum+ZP2CalibSum));
424 impPar = fbCen->Eval(ZN1CalibSum+ZP1CalibSum);
425 }
426 else if(zemenergy < (eZEMCut-deltaEZEMInf)){
427 nGenSpecNLeft = (Int_t) (fZNPer->Eval(ZN1CalibSum));
428 nGenSpecPLeft = (Int_t) (fZPPer->Eval(ZP1CalibSum));
429 nGenSpecLeft = (Int_t) (fZDCPer->Eval(ZN1CalibSum+ZP1CalibSum));
430 impPar = fbPer->Eval(ZN1CalibSum+ZP1CalibSum);
431 }
432 else if(zemenergy >= (eZEMCut-deltaEZEMInf) && zemenergy <= (eZEMCut+deltaEZEMSup)){
433 nGenSpecNLeft = (Int_t) (fZEMn->Eval(zemenergy));
434 nGenSpecPLeft = (Int_t) (fZEMp->Eval(zemenergy));
435 nGenSpecLeft = (Int_t)(fZEMsp->Eval(zemenergy));
436 impPar = fZEMb->Eval(zemenergy);
437 }
438 // ### Results from production -> 0<b<18 fm (Apr 2002)
439 if(ZN1CalibSum>162.) nGenSpecNLeft = (Int_t) (fZEMn->Eval(zemenergy));
440 if(ZP1CalibSum>59.75) nGenSpecPLeft = (Int_t) (fZEMp->Eval(zemenergy));
441 if(ZN1CalibSum+ZP1CalibSum>221.5) nGenSpecLeft = (Int_t)(fZEMsp->Eval(zemenergy));
442 if(ZN1CalibSum+ZP1CalibSum>220.) impPar = fZEMb->Eval(zemenergy);
443 */
444 //
445 //
446 // *** RECONSTRUCTION FROM REAL DATA
c0aec6f6 447 //
a4cab348 448 if(corrADCZEMHG > supValueZEM){
449 nGenSpecNLeft = (Int_t) (fZNCen->Eval(calibSumZN1HG));
450 nGenSpecPLeft = (Int_t) (fZPCen->Eval(calibSumZP1HG));
451 nGenSpecLeft = (Int_t) (fZDCCen->Eval(calibSumZN1HG+calibSumZP1HG));
452 nGenSpecNRight = (Int_t) (fZNCen->Eval(calibSumZN2HG));
453 nGenSpecPRight = (Int_t) (fZNCen->Eval(calibSumZP2HG));
454 nGenSpecRight = (Int_t) (fZNCen->Eval(calibSumZN2HG+calibSumZP2HG));
455 impPar = fbCen->Eval(calibSumZN1HG+calibSumZP1HG);
646f1679 456 }
a4cab348 457 else if(corrADCZEMHG < infValueZEM){
458 nGenSpecNLeft = (Int_t) (fZNPer->Eval(calibSumZN1HG));
459 nGenSpecPLeft = (Int_t) (fZPPer->Eval(calibSumZP1HG));
460 nGenSpecLeft = (Int_t) (fZDCPer->Eval(calibSumZN1HG+calibSumZP1HG));
461 impPar = fbPer->Eval(calibSumZN1HG+calibSumZP1HG);
646f1679 462 }
a4cab348 463 else if(corrADCZEMHG >= infValueZEM && corrADCZEMHG <= supValueZEM){
464 nGenSpecNLeft = (Int_t) (fZEMn->Eval(corrADCZEMHG));
465 nGenSpecPLeft = (Int_t) (fZEMp->Eval(corrADCZEMHG));
466 nGenSpecLeft = (Int_t)(fZEMsp->Eval(corrADCZEMHG));
467 impPar = fZEMb->Eval(corrADCZEMHG);
646f1679 468 }
469 //
a4cab348 470 if(calibSumZN1HG/maxValEZN1>1.) nGenSpecNLeft = (Int_t) (fZEMn->Eval(corrADCZEMHG));
471 if(calibSumZP1HG/maxValEZP1>1.) nGenSpecPLeft = (Int_t) (fZEMp->Eval(corrADCZEMHG));
472 if((calibSumZN1HG+calibSumZP1HG/maxValEZDC1)>1.){
473 nGenSpecLeft = (Int_t)(fZEMsp->Eval(corrADCZEMHG));
474 impPar = fZEMb->Eval(corrADCZEMHG);
646f1679 475 }
a4cab348 476 if(calibSumZN2HG/maxValEZN2>1.) nGenSpecNRight = (Int_t) (fZEMn->Eval(corrADCZEMHG));
477 if(calibSumZP2HG/maxValEZP2>1.) nGenSpecPRight = (Int_t) (fZEMp->Eval(corrADCZEMHG));
478 if((calibSumZN2HG+calibSumZP2HG/maxValEZDC2)>1.) nGenSpecRight = (Int_t)(fZEMsp->Eval(corrADCZEMHG));
646f1679 479 //
480 if(nGenSpecNLeft>125) nGenSpecNLeft=125;
481 else if(nGenSpecNLeft<0) nGenSpecNLeft=0;
482 if(nGenSpecPLeft>82) nGenSpecPLeft=82;
483 else if(nGenSpecPLeft<0) nGenSpecPLeft=0;
484 if(nGenSpecLeft>207) nGenSpecLeft=207;
485 else if(nGenSpecLeft<0) nGenSpecLeft=0;
8309c1ab 486
980685f2 487 // --- Number of generated participants (from HIJING parameterization)
646f1679 488 Int_t nPart, nPartTotLeft, nPartTotRight;
489 nPart = 207-nGenSpecNLeft-nGenSpecPLeft;
490 nPartTotLeft = 207-nGenSpecLeft;
491 nPartTotRight = 207-nGenSpecRight;
e90a5fef 492 if(nPart<0) nPart=0;
493 if(nPartTotLeft<0) nPartTotLeft=0;
494 if(nPartTotRight<0) nPartTotRight=0;
e97af564 495 //
e90a5fef 496 // *** DEBUG ***
422ed0d5 497 printf("\n\t AliZDCReconstructor -> calibSumZN1HG %1.0f, calibSumZP1HG %1.0f,"
e90a5fef 498 " calibSumZN2HG %1.0f, calibSumZP2HG %1.0f, corrADCZEMHG %1.0f\n",
499 calibSumZN1HG,calibSumZP1HG,calibSumZN2HG,calibSumZP2HG,corrADCZEMHG);
500 printf("\t AliZDCReconstructor -> nGenSpecNLeft %d, nGenSpecPLeft %d, nGenSpecLeft %d\n"
501 "\t\t nGenSpecNRight %d, nGenSpecPRight %d, nGenSpecRight %d\n",
502 nGenSpecNLeft, nGenSpecPLeft, nGenSpecLeft,
503 nGenSpecNRight, nGenSpecPRight, nGenSpecRight);
504 printf("\t AliZDCReconstructor -> NpartL %d, NpartR %d, b %1.2f fm\n\n",nPartTotLeft, nPartTotRight, impPar);
422ed0d5 505
646f1679 506 // create the output tree
a4cab348 507 AliZDCReco reco(calibSumZN1HG, calibSumZP1HG, calibSumZN2HG, calibSumZP2HG,
508 calibTowZN1LG, calibTowZN2LG, calibTowZP1LG, calibTowZP2LG,
84d6255e 509 calibTowZN1LG, calibTowZP1LG, calibTowZN2LG, calibTowZP2LG,
f5d41205 510 corrADCZEMHG,
646f1679 511 nDetSpecNLeft, nDetSpecPLeft, nDetSpecNRight, nDetSpecPRight,
512 nGenSpecNLeft, nGenSpecPLeft, nGenSpecLeft, nGenSpecNRight,
513 nGenSpecPRight, nGenSpecRight,
514 nPartTotLeft, nPartTotRight, impPar);
515
8309c1ab 516 AliZDCReco* preco = &reco;
517 const Int_t kBufferSize = 4000;
70f04f6d 518 clustersTree->Branch("ZDC", "AliZDCReco", &preco, kBufferSize);
8309c1ab 519
520 // write the output tree
70f04f6d 521 clustersTree->Fill();
8309c1ab 522}
523
524//_____________________________________________________________________________
70f04f6d 525void AliZDCReconstructor::FillZDCintoESD(TTree *clustersTree, AliESDEvent* esd) const
8309c1ab 526{
70f04f6d 527 // fill energies and number of participants to the ESD
8309c1ab 528
8309c1ab 529 AliZDCReco reco;
530 AliZDCReco* preco = &reco;
70f04f6d 531 clustersTree->SetBranchAddress("ZDC", &preco);
8309c1ab 532
70f04f6d 533 clustersTree->GetEntry(0);
84d6255e 534 //
535 esd->SetZDC(reco.GetZN1Energy(), reco.GetZP1Energy(), reco.GetZEMsignal(),
536 reco.GetZN2Energy(), reco.GetZP2Energy(),
537 reco.GetNPartLeft());
538 //
539 //
1c2caeca 540 /*Double_t tZN1Ene[4], tZN2Ene[4];
e90a5fef 541 for(Int_t i=0; i<4; i++){
542 tZN1Ene[i] = reco.GetZN1EnTow(i);
543 tZN2Ene[i] = reco.GetZN2EnTow(i);
544 }
84d6255e 545 esd->SetZN1TowerEnergy(tZN1Ene);
546 esd->SetZN2TowerEnergy(tZN2Ene);
e90a5fef 547 esd->SetZDC(tZN1Ene, tZN2Ene, reco.GetZN1Energy(), reco.GetZP1Energy(), reco.GetZEMsignal(),
646f1679 548 reco.GetZN2Energy(), reco.GetZP2Energy(),
549 reco.GetNPartLeft());
1c2caeca 550 */
a4cab348 551
8309c1ab 552}
48642b09 553
554//_____________________________________________________________________________
78d18275 555AliCDBStorage* AliZDCReconstructor::SetStorage(const char *uri)
48642b09 556{
cc2abffd 557 // Setting the storage
48642b09 558
78d18275 559 Bool_t deleteManager = kFALSE;
48642b09 560
78d18275 561 AliCDBManager *manager = AliCDBManager::Instance();
562 AliCDBStorage *defstorage = manager->GetDefaultStorage();
48642b09 563
78d18275 564 if(!defstorage || !(defstorage->Contains("ZDC"))){
565 AliWarning("No default storage set or default storage doesn't contain ZDC!");
566 manager->SetDefaultStorage(uri);
567 deleteManager = kTRUE;
568 }
569
570 AliCDBStorage *storage = manager->GetDefaultStorage();
571
572 if(deleteManager){
573 AliCDBManager::Instance()->UnsetDefaultStorage();
574 defstorage = 0; // the storage is killed by AliCDBManager::Instance()->Destroy()
575 }
576
577 return storage;
578}
48642b09 579
78d18275 580//_____________________________________________________________________________
6024ec85 581AliZDCPedestals* AliZDCReconstructor::GetPedData() const
582{
583
584 // Getting pedestal calibration object for ZDC set
585
586 AliCDBEntry *entry = AliCDBManager::Instance()->Get("ZDC/Calib/Pedestals");
587 if(!entry) AliFatal("No calibration data loaded!");
588
589 AliZDCPedestals *calibdata = dynamic_cast<AliZDCPedestals*> (entry->GetObject());
590 if(!calibdata) AliFatal("Wrong calibration object in calibration file!");
591
592 return calibdata;
593}
594
595//_____________________________________________________________________________
596AliZDCCalib* AliZDCReconstructor::GetECalibData() const
597{
598
599 // Getting energy and equalization calibration object for ZDC set
600
601 AliCDBEntry *entry = AliCDBManager::Instance()->Get("ZDC/Calib/Calib");
602 if(!entry) AliFatal("No calibration data loaded!");
603
604 AliZDCCalib *calibdata = dynamic_cast<AliZDCCalib*> (entry->GetObject());
605 if(!calibdata) AliFatal("Wrong calibration object in calibration file!");
606
607 return calibdata;
608}
609
610//_____________________________________________________________________________
611AliZDCRecParam* AliZDCReconstructor::GetRecParams() const
78d18275 612{
48642b09 613
6024ec85 614 // Getting energy and equalization calibration object for ZDC set
4fda3ba1 615
6024ec85 616 AliCDBEntry *entry = AliCDBManager::Instance()->Get("ZDC/Calib/RecParam");
457a440d 617 if(!entry) AliFatal("No calibration data loaded!");
4fda3ba1 618
6024ec85 619 AliZDCRecParam *calibdata = dynamic_cast<AliZDCRecParam*> (entry->GetObject());
457a440d 620 if(!calibdata) AliFatal("Wrong calibration object in calibration file!");
48642b09 621
78d18275 622 return calibdata;
48642b09 623}