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