]> git.uio.no Git - u/mrichter/AliRoot.git/blame_incremental - ZDC/AliZDCReconstructor.cxx
git://git.uio.no / u / mrichter / AliRoot.git / blame_incremental
? search:
summary | shortlog | log | commit | commitdiff | tree
blob | blame (non-incremental) | history | HEAD
Updated classes
[u/mrichter/AliRoot.git] / ZDC / AliZDCReconstructor.cxx
... / ...
CommitLineData
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 *
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 "AliRunInfo.h"
29#include "AliRawReader.h"
30#include "AliESDEvent.h"
31#include "AliESDZDC.h"
32#include "AliZDCDigit.h"
33#include "AliZDCRawStream.h"
34#include "AliZDCReco.h"
35#include "AliZDCReconstructor.h"
36#include "AliZDCPedestals.h"
37#include "AliZDCCalib.h"
38#include "AliZDCRecoParam.h"
39#include "AliZDCRecoParampp.h"
40#include "AliZDCRecoParamPbPb.h"
41
42
43ClassImp(AliZDCReconstructor)
44AliZDCRecoParam *AliZDCReconstructor::fRecoParam=0; //reconstruction parameters
45
46//_____________________________________________________________________________
47AliZDCReconstructor:: AliZDCReconstructor() :
48 fPedData(GetPedData()),
49 fECalibData(GetECalibData())
50{
51 // **** Default constructor
52
53}
54
55
56//_____________________________________________________________________________
57AliZDCReconstructor::~AliZDCReconstructor()
58{
59// destructor
60 if(fRecoParam) delete fRecoParam;
61 if(fPedData) delete fPedData;
62 if(fECalibData) delete fECalibData;
63}
64
65//_____________________________________________________________________________
66void AliZDCReconstructor::Reconstruct(TTree* digitsTree, TTree* clustersTree) const
67{
68 // *** Local ZDC reconstruction for digits
69 // Works on the current event
70
71 // Retrieving calibration data
72 Float_t meanPed[48];
73 for(Int_t jj=0; jj<48; jj++) meanPed[jj] = fPedData->GetMeanPed(jj);
74
75 // get digits
76 AliZDCDigit digit;
77 AliZDCDigit* pdigit = &digit;
78 digitsTree->SetBranchAddress("ZDC", &pdigit);
79 //printf("\n\t # of digits in tree: %d\n",(Int_t) digitsTree->GetEntries());
80
81 // loop over digits
82 Float_t tZN1Corr[10], tZP1Corr[10], tZN2Corr[10], tZP2Corr[10];
83 Float_t dZEM1Corr[2], dZEM2Corr[2], PMRef1[2], PMRef2[2];
84 for(Int_t i=0; i<10; i++){
85 tZN1Corr[i] = tZP1Corr[i] = tZN2Corr[i] = tZP2Corr[i] = 0.;
86 if(i<2) dZEM1Corr[i] = dZEM2Corr[i] = PMRef1[i] = PMRef2[i] = 0.;
87 }
88 //
89 for (Int_t iDigit = 0; iDigit < (digitsTree->GetEntries()/2); iDigit++) {
90 digitsTree->GetEntry(iDigit);
91 if (!pdigit) continue;
92 //
93 Int_t det = digit.GetSector(0);
94 Int_t quad = digit.GetSector(1);
95 Int_t pedindex = -1, kNch = 24;
96 //printf("\n\t Digit #%d det %d quad %d", iDigit, det, quad);
97 //
98 if(quad != 5){ // ZDC (not reference PTMs!)
99 if(det == 1){ // *** ZNC
100 pedindex = quad;
101 tZN1Corr[quad] = (Float_t) (digit.GetADCValue(0)-meanPed[pedindex]);
102 if(tZN1Corr[quad]<0.) tZN1Corr[quad] = 0.;
103 tZN1Corr[quad+5] = (Float_t) (digit.GetADCValue(1)-meanPed[pedindex+kNch]);
104 if(tZN1Corr[quad+5]<0.) tZN1Corr[quad+5] = 0.;
105 //printf("\t pedindex %d tZN1Corr[%d] = %1.0f tZN1Corr[%d] = %1.0f",
106 // pedindex, quad, tZN1Corr[quad], quad+5, tZN1Corr[quad+5]);
107 }
108 else if(det == 2){ // *** ZP1
109 pedindex = quad+5;
110 tZP1Corr[quad] = (Float_t) (digit.GetADCValue(0)-meanPed[pedindex]);
111 if(tZP1Corr[quad]<0.) tZP1Corr[quad] = 0.;
112 tZP1Corr[quad+5] = (Float_t) (digit.GetADCValue(1)-meanPed[pedindex+kNch]);
113 if(tZP1Corr[quad+5]<0.) tZP1Corr[quad+5] = 0.;
114 //printf("\t pedindex %d tZP1Corr[%d] = %1.0f tZP1Corr[%d] = %1.0f",
115 // pedindex, quad, tZP1Corr[quad], quad+5, tZP1Corr[quad+5]);
116 }
117 else if(det == 3){
118 pedindex = quad+9;
119 if(quad == 1){ // *** ZEM1
120 dZEM1Corr[0] += (Float_t) (digit.GetADCValue(0)-meanPed[pedindex]);
121 if(dZEM1Corr[0]<0.) dZEM1Corr[0] = 0.;
122 dZEM1Corr[1] += (Float_t) (digit.GetADCValue(1)-meanPed[pedindex+kNch]);
123 if(dZEM1Corr[1]<0.) dZEM1Corr[1] = 0.;
124 //printf("\t pedindex %d tZEM1Corr[%d] = %1.0f tZEM1Corr[%d] = %1.0f",
125 // pedindex, quad, tZEM1Corr[quad], quad+1, tZEM1Corr[quad+1]);
126 }
127 else if(quad == 2){ // *** ZEM2
128 dZEM2Corr[0] += (Float_t) (digit.GetADCValue(0)-meanPed[pedindex]);
129 if(dZEM2Corr[0]<0.) dZEM2Corr[0] = 0.;
130 dZEM2Corr[1] += (Float_t) (digit.GetADCValue(1)-meanPed[pedindex+kNch]);
131 if(dZEM2Corr[1]<0.) dZEM2Corr[1] = 0.;
132 //printf("\t pedindex %d tZEM2Corr[%d] = %1.0f tZEM2Corr[%d] = %1.0f",
133 // pedindex, quad, tZEM2Corr[quad], quad+1, tZEM2Corr[quad+1]);
134 }
135 }
136 else if(det == 4){ // *** ZN2
137 pedindex = quad+12;
138 tZN2Corr[quad] = (Float_t) (digit.GetADCValue(0)-meanPed[pedindex]);
139 if(tZN2Corr[quad]<0.) tZN2Corr[quad] = 0.;
140 tZN2Corr[quad+5] = (Float_t) (digit.GetADCValue(1)-meanPed[pedindex+kNch]);
141 if(tZN2Corr[quad+5]<0.) tZN2Corr[quad+5] = 0.;
142 //printf("\t pedindex %d tZN2Corr[%d] = %1.0f tZN2Corr[%d] = %1.0f\n",
143 // pedindex, quad, tZN2Corr[quad], quad+5, tZN2Corr[quad+5]);
144 }
145 else if(det == 5){ // *** ZP2
146 pedindex = quad+17;
147 tZP2Corr[quad] = (Float_t) (digit.GetADCValue(0)-meanPed[pedindex]);
148 if(tZP2Corr[quad]<0.) tZP2Corr[quad] = 0.;
149 tZP2Corr[quad+5] = (Float_t) (digit.GetADCValue(1)-meanPed[pedindex+kNch]);
150 if(tZP2Corr[quad+5]<0.) tZP2Corr[quad+5] = 0.;
151 //printf("\t pedindex %d tZP2Corr[%d] = %1.0f tZP2Corr[%d] = %1.0f\n",
152 // pedindex, quad, tZP2Corr[quad], quad+5, tZP2Corr[quad+5]);
153 }
154 }
155 else{ // Reference PMs
156 pedindex = (det-1)/3+22;
157 if(det == 1){
158 PMRef1[0] = (Float_t) (digit.GetADCValue(0)-meanPed[pedindex]);
159 if(PMRef1[0]<0.) PMRef1[0] = 0.;
160 PMRef1[1] = (Float_t) (digit.GetADCValue(1)-meanPed[pedindex+kNch]);
161 if(PMRef2[1]<0.) PMRef1[1] = 0.;
162 }
163 else if(det == 4){
164 PMRef2[0] = (Float_t) (digit.GetADCValue(0)-meanPed[pedindex]);
165 if(PMRef2[0]<0.) PMRef2[0] = 0.;
166 PMRef2[1] = (Float_t) (digit.GetADCValue(1)-meanPed[pedindex+kNch]);
167 if(PMRef2[1]<0.) PMRef2[1] = 0.;
168 }
169 }
170 }
171
172 // Reconstruct the event according to beam type from AliInfo
173 AliRunInfo *runInfo;
174 if((strcmp(runInfo->GetBeamType(),"p-p") == 0)){
175 printf("\n\t ZDC reconstructing p-p event\n");
176 ReconstructEventpp(clustersTree, tZN1Corr, tZP1Corr, tZN2Corr, tZP2Corr,
177 dZEM1Corr, dZEM2Corr, PMRef1, PMRef2);
178 }
179 else if((strcmp(runInfo->GetBeamType(),"A-A") == 0)){
180 printf("\n\t ZDC reconstructing A-A event\n");
181 ReconstructEventPbPb(clustersTree, tZN1Corr, tZP1Corr, tZN2Corr, tZP2Corr,
182 dZEM1Corr, dZEM2Corr, PMRef1, PMRef2);
183 }
184 else printf("\n\t ZDC does not reconstruct event! Beam type != p-p || A-A\n");
185}
186
187//_____________________________________________________________________________
188void AliZDCReconstructor::Reconstruct(AliRawReader* rawReader, TTree* clustersTree) const
189{
190 // *** ZDC raw data reconstruction
191 // Works on the current event
192
193 // Retrieving calibration data
194 Float_t meanPed[48];
195 for(Int_t jj=0; jj<48; jj++) meanPed[jj] = fPedData->GetMeanPed(jj);
196
197 rawReader->Reset();
198
199 // loop over raw data
200 Float_t tZN1Corr[10], tZP1Corr[10], tZN2Corr[10], tZP2Corr[10];
201 Float_t dZEM1Corr[2], dZEM2Corr[2], PMRef1[2], PMRef2[2];
202 for(Int_t i=0; i<10; i++){
203 tZN1Corr[i] = tZP1Corr[i] = tZN2Corr[i] = tZP2Corr[i] = 0.;
204 if(i<2) dZEM1Corr[i] = dZEM2Corr[i] = PMRef1[i] = PMRef2[i] = 0.;
205 }
206 //
207 AliZDCRawStream rawData(rawReader);
208 Int_t kNch = 24;
209 while (rawData.Next()) {
210 if(rawData.IsADCDataWord()){
211 Int_t det = rawData.GetSector(0);
212 Int_t quad = rawData.GetSector(1);
213 Int_t gain = rawData.GetADCGain();
214 Int_t pedindex=0;
215 //
216 if(quad !=5){ // ZDCs (not reference PTMs)
217 if(det == 1){
218 pedindex = quad;
219 if(gain == 0) tZN1Corr[quad] += (Float_t) (rawData.GetADCValue()-meanPed[pedindex]);
220 else tZN1Corr[quad+5] += (Float_t) (rawData.GetADCValue()-meanPed[pedindex+kNch]);
221 }
222 else if(det == 2){
223 pedindex = quad+5;
224 if(gain == 0) tZP1Corr[quad] += (Float_t) (rawData.GetADCValue()-meanPed[pedindex]);
225 else tZP1Corr[quad+5] += (Float_t) (rawData.GetADCValue()-meanPed[pedindex+kNch]);
226 }
227 else if(det == 3){
228 pedindex = quad+9;
229 if(quad==1){
230 if(gain == 0) dZEM1Corr[0] += (Float_t) (rawData.GetADCValue()-meanPed[pedindex]);
231 else dZEM1Corr[1] += (Float_t) (rawData.GetADCValue()-meanPed[pedindex+kNch]);
232 }
233 else if(quad==2){
234 if(gain == 0) dZEM2Corr[0] += (Float_t) (rawData.GetADCValue()-meanPed[pedindex]);
235 else dZEM2Corr[1] += (Float_t) (rawData.GetADCValue()-meanPed[pedindex+kNch]);
236 }
237 }
238 else if(det == 4){
239 pedindex = quad+12;
240 if(gain == 0) tZN2Corr[quad] += (Float_t) (rawData.GetADCValue()-meanPed[pedindex]);
241 else tZN2Corr[quad+5] += (Float_t) (rawData.GetADCValue()-meanPed[pedindex+kNch]);
242 }
243 else if(det == 5){
244 pedindex = quad+17;
245 if(gain == 0) tZP2Corr[quad] += (Float_t) (rawData.GetADCValue()-meanPed[pedindex]);
246 else tZP2Corr[quad+5] += (Float_t) (rawData.GetADCValue()-meanPed[pedindex+kNch]);
247 }
248 //printf("\t AliZDCReconstructor - det %d quad %d res %d -> Ped[%d] = %1.0f\n",
249 // det,quad,gain, pedindex, meanPed[pedindex]);
250 }
251 else{ // reference PM
252 pedindex = (det-1)/3 + 22;
253 if(det == 1){
254 if(gain==0) PMRef1[0] += (Float_t) (rawData.GetADCValue()-meanPed[pedindex]);
255 else PMRef1[1] += (Float_t) (rawData.GetADCValue()-meanPed[pedindex]);
256 }
257 else if(det ==4){
258 if(gain==0) PMRef2[0] += (Float_t) (rawData.GetADCValue()-meanPed[pedindex]);
259 else PMRef2[1] += (Float_t) (rawData.GetADCValue()-meanPed[pedindex]);
260 }
261 }
262 }//IsADCDataWord
263 }
264
265 // Reconstruct the event according to beam type from AliInfo
266 AliRunInfo *runInfo;
267 if((strcmp(runInfo->GetBeamType(),"p-p") == 0)){
268 printf("\n\t ZDC reconstructing p-p event\n");
269 ReconstructEventpp(clustersTree, tZN1Corr, tZP1Corr, tZN2Corr, tZP2Corr,
270 dZEM1Corr, dZEM2Corr, PMRef1, PMRef2);
271 }
272 else if((strcmp(runInfo->GetBeamType(),"A-A") == 0)){
273 printf("\n\t ZDC reconstructing A-A event\n");
274 ReconstructEventPbPb(clustersTree, tZN1Corr, tZP1Corr, tZN2Corr, tZP2Corr,
275 dZEM1Corr, dZEM2Corr, PMRef1, PMRef2);
276 }
277 else printf("\n\t ZDC does not reconstruct event! Beam type != p-p || A-A\n");
278
279}
280
281//_____________________________________________________________________________
282void AliZDCReconstructor::ReconstructEventpp(TTree *clustersTree, Float_t* ZN1ADCCorr,
283 Float_t* ZP1ADCCorr, Float_t* ZN2ADCCorr, Float_t* ZP2ADCCorr,
284 Float_t* ZEM1ADCCorr, Float_t* ZEM2ADCCorr, Float_t* PMRef1, Float_t* PMRef2) const
285{
286 // ***** Reconstruct one event
287
288 // *** RECONSTRUCTION FROM "REAL" DATA
289 //
290 // Retrieving calibration data
291 // --- Equalization coefficients ---------------------------------------------
292 Float_t equalCoeffZN1[5], equalCoeffZP1[5], equalCoeffZN2[5], equalCoeffZP2[5];
293 for(Int_t ji=0; ji<5; ji++){
294 equalCoeffZN1[ji] = fECalibData->GetZN1EqualCoeff(ji);
295 equalCoeffZP1[ji] = fECalibData->GetZP1EqualCoeff(ji);
296 equalCoeffZN2[ji] = fECalibData->GetZN2EqualCoeff(ji);
297 equalCoeffZP2[ji] = fECalibData->GetZP2EqualCoeff(ji);
298 }
299 // --- Energy calibration factors ------------------------------------
300 Float_t calibEne[4];
301 // *********************************************************************
302 // **** Until the beam type info isn't known @ reconstruction level ****
303 // **** the energy calibration coefficient are manually set to 1 ****
304 // **** as it will be in real life for pp data taking ****
305 // *********************************************************************
306 //for(Int_t ij=0; ij<4; ij++) calibEne[ij] = fECalibData->GetEnCalib(ij);
307 for(Int_t ij=0; ij<4; ij++) calibEne[ij] = 1.;
308
309 // Equalization of detector responses
310 Float_t equalTowZN1[10], equalTowZN2[10], equalTowZP1[10], equalTowZP2[10];
311 for(Int_t gi=0; gi<5; gi++){
312 equalTowZN1[gi] = ZN1ADCCorr[gi]*equalCoeffZN1[gi];
313 equalTowZN1[gi+5] = ZN1ADCCorr[gi+5]*equalCoeffZN1[gi];
314 equalTowZP1[gi] = ZP1ADCCorr[gi]*equalCoeffZP1[gi];
315 equalTowZP1[gi+5] = ZP1ADCCorr[gi+5]*equalCoeffZP1[gi];
316 equalTowZN2[gi] = ZN2ADCCorr[gi]*equalCoeffZN2[gi];
317 equalTowZN2[gi+5] = ZN2ADCCorr[gi+5]*equalCoeffZN2[gi];
318 equalTowZP2[gi] = ZP2ADCCorr[gi]*equalCoeffZP2[gi];
319 equalTowZP2[gi+5] = ZP2ADCCorr[gi+5]*equalCoeffZP2[gi];
320 }
321
322 // Energy calibration of detector responses
323 Float_t calibTowZN1[10], calibTowZN2[10], calibTowZP1[10], calibTowZP2[10];
324 Float_t calibSumZN1[]={0,0}, calibSumZN2[]={0,0}, calibSumZP1[]={0,0}, calibSumZP2[]={0,0};
325 for(Int_t gi=0; gi<10; gi++){
326 calibTowZN1[gi] = equalTowZN1[gi]*calibEne[0];
327 calibTowZP1[gi] = equalTowZP1[gi]*calibEne[1];
328 calibTowZN2[gi] = equalTowZN2[gi]*calibEne[2];
329 calibTowZP2[gi] = equalTowZP2[gi]*calibEne[3];
330 //
331 if(gi<5){
332 calibSumZN1[0] += calibTowZN1[gi];
333 calibSumZP1[0] += calibTowZP1[gi];
334 calibSumZN2[0] += calibTowZN2[gi];
335 calibSumZP2[0] += calibTowZP2[gi];
336 }
337 //
338 else{
339 calibSumZN1[1] += calibTowZN1[gi];
340 calibSumZP1[1] += calibTowZP1[gi];
341 calibSumZN2[1] += calibTowZN2[gi];
342 calibSumZP2[1] += calibTowZP2[gi];
343 }
344 }
345
346 //
347 // --- Reconstruction parameters ------------------
348 if(!fRecoParam) fRecoParam = (AliZDCRecoParampp*) AliZDCRecoParampp::GetppRecoParam();
349
350 // --- Number of detected spectator nucleons
351 // *** N.B. -> It works only in Pb-Pb
352 Int_t nDetSpecNLeft, nDetSpecPLeft, nDetSpecNRight, nDetSpecPRight;
353 Float_t beamE = fRecoParam->GetBeamEnergy();
354 nDetSpecNLeft = (Int_t) (calibSumZN1[0]/beamE);
355 nDetSpecPLeft = (Int_t) (calibSumZP1[0]/beamE);
356 nDetSpecNRight = (Int_t) (calibSumZN2[0]/beamE);
357 nDetSpecPRight = (Int_t) (calibSumZP2[0]/beamE);
358 /*printf("\n\t AliZDCReconstructor -> nDetSpecNLeft %d, nDetSpecPLeft %d,"
359 " nDetSpecNRight %d, nDetSpecPRight %d\n",nDetSpecNLeft, nDetSpecPLeft,
360 nDetSpecNRight, nDetSpecPRight);*/
361
362 // --- Number of generated spectator nucleons (from HIJING parameterization)
363 Int_t nGenSpecNLeft=0, nGenSpecPLeft=0, nGenSpecLeft=0;
364 Int_t nGenSpecNRight=0, nGenSpecPRight=0, nGenSpecRight=0;
365 Int_t nPartTotLeft=0, nPartTotRight=0;
366 Double_t impPar=0.;
367
368 // create the output tree
369 AliZDCReco reco(calibSumZN1, calibSumZP1, calibSumZN2, calibSumZP2,
370 calibTowZN1, calibTowZP1, calibTowZN2, calibTowZP2,
371 ZEM1ADCCorr, ZEM2ADCCorr, PMRef1, PMRef2,
372 nDetSpecNLeft, nDetSpecPLeft, nDetSpecNRight, nDetSpecPRight,
373 nGenSpecNLeft, nGenSpecPLeft, nGenSpecLeft, nGenSpecNRight,
374 nGenSpecPRight, nGenSpecRight, nPartTotLeft, nPartTotRight, impPar);
375
376 AliZDCReco* preco = &reco;
377 const Int_t kBufferSize = 4000;
378 clustersTree->Branch("ZDC", "AliZDCReco", &preco, kBufferSize);
379
380 // write the output tree
381 clustersTree->Fill();
382}
383
384//_____________________________________________________________________________
385void AliZDCReconstructor::ReconstructEventPbPb(TTree *clustersTree, Float_t* ZN1ADCCorr,
386 Float_t* ZP1ADCCorr, Float_t* ZN2ADCCorr, Float_t* ZP2ADCCorr,
387 Float_t* ZEM1ADCCorr, Float_t* ZEM2ADCCorr, Float_t* PMRef1, Float_t* PMRef2) const
388{
389 // ***** Reconstruct one event
390
391 // *** RECONSTRUCTION FROM "REAL" DATA
392 //
393 // Retrieving calibration data
394 // --- Equalization coefficients ---------------------------------------------
395 Float_t equalCoeffZN1[5], equalCoeffZP1[5], equalCoeffZN2[5], equalCoeffZP2[5];
396 for(Int_t ji=0; ji<5; ji++){
397 equalCoeffZN1[ji] = fECalibData->GetZN1EqualCoeff(ji);
398 equalCoeffZP1[ji] = fECalibData->GetZP1EqualCoeff(ji);
399 equalCoeffZN2[ji] = fECalibData->GetZN2EqualCoeff(ji);
400 equalCoeffZP2[ji] = fECalibData->GetZP2EqualCoeff(ji);
401 }
402 // --- Energy calibration factors ------------------------------------
403 Float_t calibEne[4];
404 for(Int_t ij=0; ij<4; ij++) calibEne[ij] = fECalibData->GetEnCalib(ij);
405
406 // Equalization of detector responses
407 Float_t equalTowZN1[10], equalTowZN2[10], equalTowZP1[10], equalTowZP2[10];
408 for(Int_t gi=0; gi<5; gi++){
409 equalTowZN1[gi] = ZN1ADCCorr[gi]*equalCoeffZN1[gi];
410 equalTowZN1[gi+5] = ZN1ADCCorr[gi+5]*equalCoeffZN1[gi];
411 equalTowZP1[gi] = ZP1ADCCorr[gi]*equalCoeffZP1[gi];
412 equalTowZP1[gi+5] = ZP1ADCCorr[gi+5]*equalCoeffZP1[gi];
413 equalTowZN2[gi] = ZN2ADCCorr[gi]*equalCoeffZN2[gi];
414 equalTowZN2[gi+5] = ZN2ADCCorr[gi+5]*equalCoeffZN2[gi];
415 equalTowZP2[gi] = ZP2ADCCorr[gi]*equalCoeffZP2[gi];
416 equalTowZP2[gi+5] = ZP2ADCCorr[gi+5]*equalCoeffZP2[gi];
417 }
418
419 // Energy calibration of detector responses
420 Float_t calibTowZN1[10], calibTowZN2[10], calibTowZP1[10], calibTowZP2[10];
421 Float_t calibSumZN1[]={0,0}, calibSumZN2[]={0,0}, calibSumZP1[]={0,0}, calibSumZP2[]={0,0};
422 for(Int_t gi=0; gi<10; gi++){
423 calibTowZN1[gi] = equalTowZN1[gi]*calibEne[0];
424 calibTowZP1[gi] = equalTowZP1[gi]*calibEne[1];
425 calibTowZN2[gi] = equalTowZN2[gi]*calibEne[2];
426 calibTowZP2[gi] = equalTowZP2[gi]*calibEne[3];
427 //
428 if(gi<5){
429 calibSumZN1[0] += calibTowZN1[gi];
430 calibSumZP1[0] += calibTowZP1[gi];
431 calibSumZN2[0] += calibTowZN2[gi];
432 calibSumZP2[0] += calibTowZP2[gi];
433 }
434 //
435 else{
436 calibSumZN1[1] += calibTowZN1[gi];
437 calibSumZP1[1] += calibTowZP1[gi];
438 calibSumZN2[1] += calibTowZN2[gi];
439 calibSumZP2[1] += calibTowZP2[gi];
440 }
441 }
442
443 //
444 // --- Reconstruction parameters ------------------
445 if(!fRecoParam) fRecoParam = (AliZDCRecoParamPbPb*) AliZDCRecoParamPbPb::GetPbPbRecoParam();
446 //
447 Float_t endPointZEM = fRecoParam->GetZEMEndValue();
448 Float_t cutFractionZEM = fRecoParam->GetZEMCutFraction();
449 Float_t dZEMSup = fRecoParam->GetDZEMSup();
450 Float_t dZEMInf = fRecoParam->GetDZEMInf();
451 //
452 Float_t cutValueZEM = endPointZEM*cutFractionZEM;
453 Float_t supValueZEM = cutValueZEM+(endPointZEM*dZEMSup);
454 Float_t infValueZEM = cutValueZEM-(endPointZEM*dZEMInf);
455 //
456 Float_t maxValEZN1 = fRecoParam->GetEZN1MaxValue();
457 Float_t maxValEZP1 = fRecoParam->GetEZP1MaxValue();
458 Float_t maxValEZDC1 = fRecoParam->GetEZDC1MaxValue();
459 Float_t maxValEZN2 = fRecoParam->GetEZN2MaxValue();
460 Float_t maxValEZP2 = fRecoParam->GetEZP2MaxValue();
461 Float_t maxValEZDC2 = fRecoParam->GetEZDC2MaxValue();
462 //
463 //printf("\n\t AliZDCReconstructor -> ZEMEndPoint %1.0f, ZEMCutValue %1.0f,"
464 // " ZEMSupValue %1.0f, ZEMInfValue %1.0f\n",endPointZEM,cutValueZEM,supValueZEM,infValueZEM);
465
466 // --- Number of detected spectator nucleons
467 // *** N.B. -> It works only in Pb-Pb
468 Int_t nDetSpecNLeft, nDetSpecPLeft, nDetSpecNRight, nDetSpecPRight;
469 Float_t beamE = fRecoParam->GetBeamEnergy();
470 nDetSpecNLeft = (Int_t) (calibSumZN1[0]/beamE);
471 nDetSpecPLeft = (Int_t) (calibSumZP1[0]/beamE);
472 nDetSpecNRight = (Int_t) (calibSumZN2[0]/beamE);
473 nDetSpecPRight = (Int_t) (calibSumZP2[0]/beamE);
474 /*printf("\n\t AliZDCReconstructor -> nDetSpecNLeft %d, nDetSpecPLeft %d,"
475 " nDetSpecNRight %d, nDetSpecPRight %d\n",nDetSpecNLeft, nDetSpecPLeft,
476 nDetSpecNRight, nDetSpecPRight);*/
477
478 // --- Number of generated spectator nucleons (from HIJING parameterization)
479 Int_t nGenSpecNLeft=0, nGenSpecPLeft=0, nGenSpecLeft=0;
480 Int_t nGenSpecNRight=0, nGenSpecPRight=0, nGenSpecRight=0;
481 Double_t impPar=0.;
482 //
483 Float_t corrADCZEMHG = ZEM1ADCCorr[0] + ZEM2ADCCorr[0];
484 //
485 if(corrADCZEMHG > supValueZEM){
486 nGenSpecNLeft = (Int_t) ((fRecoParam->GetfZNCen())->Eval(calibSumZN1[0]));
487 nGenSpecPLeft = (Int_t) ((fRecoParam->GetfZPCen())->Eval(calibSumZP1[0]));
488 nGenSpecLeft = (Int_t) ((fRecoParam->GetfZDCCen())->Eval(calibSumZN1[0]+calibSumZP1[0]));
489 nGenSpecNRight = (Int_t) ((fRecoParam->GetfZNCen())->Eval(calibSumZN2[0]));
490 nGenSpecPRight = (Int_t) ((fRecoParam->GetfZNCen())->Eval(calibSumZP2[0]));
491 nGenSpecRight = (Int_t) ((fRecoParam->GetfZNCen())->Eval(calibSumZN2[0]+calibSumZP2[0]));
492 impPar = (fRecoParam->GetfbCen())->Eval(calibSumZN1[0]+calibSumZP1[0]);
493 }
494 else if(corrADCZEMHG < infValueZEM){
495 nGenSpecNLeft = (Int_t) ((fRecoParam->GetfZNPer())->Eval(calibSumZN1[0]));
496 nGenSpecPLeft = (Int_t) ((fRecoParam->GetfZPPer())->Eval(calibSumZP1[0]));
497 nGenSpecLeft = (Int_t) ((fRecoParam->GetfZDCPer())->Eval(calibSumZN1[0]+calibSumZP1[0]));
498 impPar = (fRecoParam->GetfbPer())->Eval(calibSumZN1[0]+calibSumZP1[0]);
499 }
500 else if(corrADCZEMHG >= infValueZEM && corrADCZEMHG <= supValueZEM){
501 nGenSpecNLeft = (Int_t) ((fRecoParam->GetfZEMn())->Eval(corrADCZEMHG));
502 nGenSpecPLeft = (Int_t) ((fRecoParam->GetfZEMp())->Eval(corrADCZEMHG));
503 nGenSpecLeft = (Int_t)((fRecoParam->GetfZEMsp())->Eval(corrADCZEMHG));
504 impPar = (fRecoParam->GetfZEMb())->Eval(corrADCZEMHG);
505 }
506 //
507 if(calibSumZN1[0]/maxValEZN1>1.) nGenSpecNLeft = (Int_t) ((fRecoParam->GetfZEMn())->Eval(corrADCZEMHG));
508 if(calibSumZP1[0]/maxValEZP1>1.) nGenSpecPLeft = (Int_t) ((fRecoParam->GetfZEMp())->Eval(corrADCZEMHG));
509 if((calibSumZN1[0]+calibSumZP1[0]/maxValEZDC1)>1.){
510 nGenSpecLeft = (Int_t)((fRecoParam->GetfZEMsp())->Eval(corrADCZEMHG));
511 impPar = (fRecoParam->GetfZEMb())->Eval(corrADCZEMHG);
512 }
513 if(calibSumZN2[0]/maxValEZN2>1.) nGenSpecNRight = (Int_t) ((fRecoParam->GetfZEMn())->Eval(corrADCZEMHG));
514 if(calibSumZP2[0]/maxValEZP2>1.) nGenSpecPRight = (Int_t) ((fRecoParam->GetfZEMp())->Eval(corrADCZEMHG));
515 if((calibSumZN2[0]+calibSumZP2[0]/maxValEZDC2)>1.) nGenSpecRight = (Int_t)((fRecoParam->GetfZEMsp())->Eval(corrADCZEMHG));
516 //
517 if(nGenSpecNLeft>125) nGenSpecNLeft=125;
518 else if(nGenSpecNLeft<0) nGenSpecNLeft=0;
519 if(nGenSpecPLeft>82) nGenSpecPLeft=82;
520 else if(nGenSpecPLeft<0) nGenSpecPLeft=0;
521 if(nGenSpecLeft>207) nGenSpecLeft=207;
522 else if(nGenSpecLeft<0) nGenSpecLeft=0;
523
524 // --- Number of generated participants (from HIJING parameterization)
525 Int_t nPart, nPartTotLeft, nPartTotRight;
526 nPart = 207-nGenSpecNLeft-nGenSpecPLeft;
527 nPartTotLeft = 207-nGenSpecLeft;
528 nPartTotRight = 207-nGenSpecRight;
529 if(nPart<0) nPart=0;
530 if(nPartTotLeft<0) nPartTotLeft=0;
531 if(nPartTotRight<0) nPartTotRight=0;
532 //
533 // *** DEBUG ***
534 /*printf("\n\t AliZDCReconstructor -> calibSumZN1[0] %1.0f, calibSumZP1[0] %1.0f,"
535 " calibSumZN2[0] %1.0f, calibSumZP2[0] %1.0f, corrADCZEMHG %1.0f\n",
536 calibSumZN1[0],calibSumZP1[0],calibSumZN2[0],calibSumZP2[0],corrADCZEMHG);
537 printf("\t AliZDCReconstructor -> nGenSpecNLeft %d, nGenSpecPLeft %d, nGenSpecLeft %d\n"
538 "\t\t nGenSpecNRight %d, nGenSpecPRight %d, nGenSpecRight %d\n",
539 nGenSpecNLeft, nGenSpecPLeft, nGenSpecLeft,
540 nGenSpecNRight, nGenSpecPRight, nGenSpecRight);
541 printf("\t AliZDCReconstructor -> NpartL %d, NpartR %d, b %1.2f fm\n\n",nPartTotLeft, nPartTotRight, impPar);
542 */
543
544 // create the output tree
545 AliZDCReco reco(calibSumZN1, calibSumZP1, calibSumZN2, calibSumZP2,
546 calibTowZN1, calibTowZP1, calibTowZN2, calibTowZP2,
547 ZEM1ADCCorr, ZEM2ADCCorr, PMRef1, PMRef2,
548 nDetSpecNLeft, nDetSpecPLeft, nDetSpecNRight, nDetSpecPRight,
549 nGenSpecNLeft, nGenSpecPLeft, nGenSpecLeft, nGenSpecNRight,
550 nGenSpecPRight, nGenSpecRight, nPartTotLeft, nPartTotRight, impPar);
551
552 AliZDCReco* preco = &reco;
553 const Int_t kBufferSize = 4000;
554 clustersTree->Branch("ZDC", "AliZDCReco", &preco, kBufferSize);
555
556 // write the output tree
557 clustersTree->Fill();
558}
559
560//_____________________________________________________________________________
561void AliZDCReconstructor::FillZDCintoESD(TTree *clustersTree, AliESDEvent* esd) const
562{
563 // fill energies and number of participants to the ESD
564
565 AliZDCReco reco;
566 AliZDCReco* preco = &reco;
567 clustersTree->SetBranchAddress("ZDC", &preco);
568
569 clustersTree->GetEntry(0);
570 //
571 AliESDZDC * esdzdc = esd->GetESDZDC();
572 Float_t tZN1Ene[5], tZN2Ene[5], tZP1Ene[5], tZP2Ene[5];
573 Float_t tZN1EneLR[5], tZN2EneLR[5], tZP1EneLR[5], tZP2EneLR[5];
574 for(Int_t i=0; i<5; i++){
575 tZN1Ene[i] = reco.GetZN1HREnTow(i);
576 tZN2Ene[i] = reco.GetZN2HREnTow(i);
577 tZP1Ene[i] = reco.GetZP1HREnTow(i);
578 tZP2Ene[i] = reco.GetZP2HREnTow(i);
579 //
580 tZN1EneLR[i] = reco.GetZN1LREnTow(i);
581 tZN2EneLR[i] = reco.GetZN2LREnTow(i);
582 tZP1EneLR[i] = reco.GetZP1LREnTow(i);
583 tZP2EneLR[i] = reco.GetZP2LREnTow(i);
584 }
585 esdzdc->SetZN1TowerEnergy(tZN1Ene);
586 esdzdc->SetZN2TowerEnergy(tZN2Ene);
587 esdzdc->SetZP1TowerEnergy(tZP1Ene);
588 esdzdc->SetZP2TowerEnergy(tZP2Ene);
589 esdzdc->SetZN1TowerEnergyLR(tZN1EneLR);
590 esdzdc->SetZN2TowerEnergyLR(tZN2EneLR);
591 esdzdc->SetZP1TowerEnergyLR(tZP1EneLR);
592 esdzdc->SetZP2TowerEnergyLR(tZP2EneLR);
593 //
594 esd->SetZDC(reco.GetZN1HREnergy(), reco.GetZP1HREnergy(), reco.GetZEM1HRsignal(),
595 reco.GetZEM2HRsignal(), reco.GetZN2HREnergy(), reco.GetZP2HREnergy(),
596 reco.GetNPartLeft(),reco.GetNPartRight());
597 //
598
599}
600
601//_____________________________________________________________________________
602AliCDBStorage* AliZDCReconstructor::SetStorage(const char *uri)
603{
604 // Setting the storage
605
606 Bool_t deleteManager = kFALSE;
607
608 AliCDBManager *manager = AliCDBManager::Instance();
609 AliCDBStorage *defstorage = manager->GetDefaultStorage();
610
611 if(!defstorage || !(defstorage->Contains("ZDC"))){
612 AliWarning("No default storage set or default storage doesn't contain ZDC!");
613 manager->SetDefaultStorage(uri);
614 deleteManager = kTRUE;
615 }
616
617 AliCDBStorage *storage = manager->GetDefaultStorage();
618
619 if(deleteManager){
620 AliCDBManager::Instance()->UnsetDefaultStorage();
621 defstorage = 0; // the storage is killed by AliCDBManager::Instance()->Destroy()
622 }
623
624 return storage;
625}
626
627//_____________________________________________________________________________
628AliZDCPedestals* AliZDCReconstructor::GetPedData() const
629{
630
631 // Getting pedestal calibration object for ZDC set
632
633 AliCDBEntry *entry = AliCDBManager::Instance()->Get("ZDC/Calib/Pedestals");
634 if(!entry) AliFatal("No calibration data loaded!");
635
636 AliZDCPedestals *calibdata = dynamic_cast<AliZDCPedestals*> (entry->GetObject());
637 if(!calibdata) AliFatal("Wrong calibration object in calibration file!");
638
639 return calibdata;
640}
641
642//_____________________________________________________________________________
643AliZDCCalib* AliZDCReconstructor::GetECalibData() const
644{
645
646 // Getting energy and equalization calibration object for ZDC set
647
648 AliCDBEntry *entry = AliCDBManager::Instance()->Get("ZDC/Calib/EMDCalib");
649 if(!entry) AliFatal("No calibration data loaded!");
650
651 AliZDCCalib *calibdata = dynamic_cast<AliZDCCalib*> (entry->GetObject());
652 if(!calibdata) AliFatal("Wrong calibration object in calibration file!");
653
654 return calibdata;
655}
656