]>
Commit | Line | Data |
---|---|---|
65abd48b | 1 | #include "AliDisplacedVertexSelection.h" |
2 | #include <iostream> | |
3 | #include <TROOT.h> | |
4 | #include "AliESDEvent.h" | |
5 | #include "AliESDZDC.h" | |
6 | ClassImp(AliDisplacedVertexSelection) | |
7 | #if 0 | |
8 | ; // This is for Emacs | |
9 | #endif | |
10 | ||
11 | //____________________________________________________________________ | |
12 | AliDisplacedVertexSelection::AliDisplacedVertexSelection() | |
241cca4d | 13 | : TObject(), |
14 | fVertexZ(9999), | |
15 | fCent(100) | |
65abd48b | 16 | { |
17 | } | |
18 | //____________________________________________________________________ | |
19 | AliDisplacedVertexSelection::AliDisplacedVertexSelection(const AliDisplacedVertexSelection& o) | |
241cca4d | 20 | : TObject(o), |
21 | fVertexZ(9999), | |
22 | fCent(100) | |
65abd48b | 23 | { |
24 | } | |
25 | //____________________________________________________________________ | |
26 | AliDisplacedVertexSelection& | |
27 | AliDisplacedVertexSelection::operator=(const AliDisplacedVertexSelection& o) | |
28 | { | |
29 | if (&o == this) return *this; | |
30 | return *this; | |
31 | } | |
32 | ||
33 | //____________________________________________________________________ | |
34 | void | |
35 | AliDisplacedVertexSelection::Output(TList* /*l*/, const char* /* name*/) const | |
36 | { | |
37 | } | |
38 | ||
39 | //____________________________________________________________________ | |
40 | void | |
41 | AliDisplacedVertexSelection::Print(Option_t*) const | |
42 | { | |
6ab100ec | 43 | #if 0 |
65abd48b | 44 | char ind[gROOT->GetDirLevel()+1]; |
45 | for (Int_t i = 0; i < gROOT->GetDirLevel(); i++) ind[i] = ' '; | |
46 | ind[gROOT->GetDirLevel()] = '\0'; | |
47 | std::cout << std::boolalpha | |
48 | << std::noboolalpha << std::endl; | |
6ab100ec | 49 | #endif |
65abd48b | 50 | } |
241cca4d | 51 | |
65abd48b | 52 | //____________________________________________________________________ |
241cca4d | 53 | Bool_t |
54 | AliDisplacedVertexSelection::Process(const AliESDEvent* esd) | |
55 | { | |
56 | fVertexZ = 9999; // Default vertex value | |
57 | fCent = 100; // Default centrality value | |
58 | ||
59 | // Some constants | |
60 | const Float_t kZDCrefSum = -66.5; | |
61 | const Float_t kZDCrefDelta = -2.10; | |
62 | const Float_t kZDCsigmaSum = 3.25; | |
63 | const Float_t kZDCsigmaDelta = 2.25; | |
64 | const Float_t kZDCsigmaSumSat = 2.50; | |
65 | const Float_t kZDCsigmaDeltaSat = 1.35; | |
66 | ||
67 | // Corrections for magnetic fields | |
68 | const Float_t kZEMcorrPlusPlus[21] = {0.6840,0.7879,0.8722, | |
69 | 0.9370,0.9837,1.0137, | |
70 | 1.0292,1.0327,1.0271, | |
71 | 1.0152,1.0000,0.9844, | |
72 | 0.9714,0.9634,0.9626, | |
73 | 0.9708,0.9891,1.0181, | |
74 | 1.0574,1.1060,1.1617}; | |
75 | const Float_t kZEMcorrMoinsMoins[21] = {0.7082,0.8012,0.8809, | |
76 | 0.9447,0.9916,1.0220, | |
77 | 1.0372,1.0395,1.0318, | |
78 | 1.0174,1.0000,0.9830, | |
79 | 0.9699,0.9635,0.9662, | |
80 | 0.9794,1.0034,1.0371, | |
81 | 1.0782,1.1224,1.1634}; | |
82 | ||
83 | // --- Get the ESD object ------------------------------------------ | |
84 | AliESDZDC* esdZDC = esd->GetESDZDC(); | |
85 | if (!esdZDC) { | |
86 | AliWarning("No ZDC ESD object!"); | |
87 | return false; | |
88 | } | |
89 | Double_t currentL3 = esd->GetCurrentL3(); | |
90 | Double_t currentDipo = esd->GetCurrentDip(); | |
91 | ||
92 | // --- ZDC and ZEM energy signal ----------------------------------- | |
93 | Double_t zdcEn = (esdZDC->GetZDCN1Energy()+ | |
94 | esdZDC->GetZDCP1Energy()+ | |
95 | esdZDC->GetZDCN2Energy()+ | |
96 | esdZDC->GetZDCP2Energy()); | |
97 | Double_t zemEn = (esdZDC->GetZDCEMEnergy(0)+ | |
98 | esdZDC->GetZDCEMEnergy(1))/8.; | |
99 | ||
100 | // --- Calculate DeltaT and sumT ----------------------------------- | |
101 | Double_t deltaTdc = 0; | |
102 | Double_t sumTdc = 0; | |
103 | ||
104 | for(Int_t i = 0; i < 4; ++i) { | |
105 | if(esdZDC->GetZDCTDCData(10,i) != 0) { | |
106 | Double_t tdcCnoCorr = 0.025*(esdZDC->GetZDCTDCData(10,i)- | |
107 | esdZDC->GetZDCTDCData(14,i)); | |
108 | Double_t tdcC = esdZDC->GetZDCTDCCorrected(10,i); | |
109 | for(Int_t j = 0; j < 4; ++j) { | |
110 | if(esdZDC->GetZDCTDCData(12,j) != 0) { | |
111 | Double_t tdcAnoCorr = 0.025*(esdZDC->GetZDCTDCData(12,j)- | |
112 | esdZDC->GetZDCTDCData(14,j)); | |
113 | Double_t tdcA = esdZDC->GetZDCTDCCorrected(12,j); | |
114 | if(esdZDC->TestBit(AliESDZDC::kCorrectedTDCFilled)) { | |
115 | deltaTdc = tdcC-tdcA; | |
116 | sumTdc = tdcC+tdcA; | |
117 | } | |
118 | else { | |
119 | deltaTdc = tdcCnoCorr-tdcAnoCorr; | |
120 | sumTdc = tdcCnoCorr+tdcAnoCorr; | |
121 | } | |
122 | } | |
123 | } | |
124 | } | |
125 | } | |
126 | ||
127 | // --- Find the vertex --------------------------------------------- | |
128 | if(deltaTdc!=0. || sumTdc!=0.) { | |
129 | for (Int_t k = -10; k <= 10; ++k) { | |
130 | Float_t zsat = 2.55F * k; | |
131 | Float_t delta = (k == 0 ? kZDCsigmaDelta : kZDCsigmaDeltaSat); | |
132 | Float_t sum = (k == 0 ? kZDCsigmaSum : kZDCsigmaSumSat); | |
133 | Float_t dT = deltaTdc - kZDCrefDelta - zsat; | |
134 | Float_t sT = sumTdc - kZDCrefSum - zsat; | |
135 | Float_t check = dT * dT / delta / delta + sT * sT / sum / sum; | |
136 | if (check > 1.0 || k == 0) continue; | |
137 | ||
138 | // Set the vertex | |
139 | fVertexZ = 37.5 * k; | |
140 | ||
141 | // Correct zem energy | |
142 | if(currentDipo>0 && currentL3>0) zemEn /= kZEMcorrPlusPlus[k+10]; | |
143 | if(currentDipo<0 && currentL3<0) zemEn /= kZEMcorrMoinsMoins[k+10]; | |
144 | } | |
145 | } | |
146 | ||
147 | // --- Calculate the centrality ------------------------------------ | |
148 | Float_t c1, c2, c3, c4, c5; | |
149 | Int_t runnumber = esd->GetRunNumber(); | |
150 | if (runnumber < 137722 && runnumber >= 137161 ) { | |
151 | c1 = 16992; | |
152 | c2 = 345; | |
153 | c3 = 2.23902e+02; | |
154 | c4 = 1.56667; | |
155 | c5 = 9.49434e-05; | |
156 | } | |
157 | else if (runnumber < 139172 && runnumber >= 137722) { | |
158 | c1 = 15000; | |
159 | c2 = 295; | |
160 | c3 = 2.23660e+02; | |
161 | c4 = 1.56664; | |
162 | c5 = 8.99571e-05; | |
163 | } | |
164 | else if (runnumber >= 139172) { | |
165 | c1 = 16992; | |
166 | c2 = 345; | |
167 | c3 = 2.04789e+02; | |
168 | c4 = 1.56629; | |
169 | c5 = 1.02768e-04; | |
170 | } | |
171 | else { | |
172 | c1 = 16992; | |
173 | c2 = 345; | |
174 | c3 = 2.04789e+02; | |
175 | c4 = 1.56629; | |
176 | c5 = 1.02768e-04; | |
177 | } | |
178 | ||
179 | if (zemEn > c2) { | |
180 | Float_t slope = (zdcEn + c1) / (zemEn - c2) + c3; | |
181 | Float_t zdcCent = (TMath::ATan(slope) - c4) / c5; | |
182 | if (zdcCent >= 0) fCent = zdcCent; | |
183 | } | |
65abd48b | 184 | |
241cca4d | 185 | return true; |
186 | } | |
187 | ||
188 | //____________________________________________________________________ | |
65abd48b | 189 | Double_t |
190 | AliDisplacedVertexSelection::CheckDisplacedVertex(const AliESDEvent* esd) const { | |
191 | // Code taken directly from M.Guilbaud. | |
241cca4d | 192 | AliWarning("This method is deprecated"); |
65abd48b | 193 | Double_t zvtx = 9999.; |
194 | ||
195 | Float_t kZDCrefSum =-66.5; | |
196 | Float_t kZDCrefDelta =-2.10; | |
197 | Float_t kZDCsigmaSum = 3.25; | |
198 | Float_t kZDCsigmaDelta = 2.25; | |
199 | Float_t kZDCsigmaSumSat = 2.50; | |
200 | Float_t kZDCsigmaDeltaSat = 1.35; | |
201 | ||
202 | AliESDZDC* esdZDC = esd -> GetESDZDC(); | |
203 | ||
204 | /* Double_t zdcNCEnergy = esdZDC->GetZDCN1Energy(); | |
205 | Double_t fZpcEnergy = esdZDC->GetZDCP1Energy(); | |
206 | Double_t fZnaEnergy = esdZDC->GetZDCN2Energy(); | |
207 | Double_t fZpaEnergy = esdZDC->GetZDCP2Energy(); | |
208 | Double_t fZem1Energy = esdZDC->GetZDCEMEnergy(0); | |
209 | Double_t fZem2Energy = esdZDC->GetZDCEMEnergy(1);*/ | |
210 | ||
211 | //Double_t fzdcEn = (esdZDC->GetZDCN1Energy()+esdZDC->GetZDCP1Energy()+esdZDC->GetZDCN2Energy()+esdZDC->GetZDCP2Energy()); | |
212 | //Double_t fzemEn = (esdZDC->GetZDCEMEnergy(0)+esdZDC->GetZDCEMEnergy(1))/8.; | |
213 | ||
214 | ||
241cca4d | 215 | |
216 | /////////////////// | |
217 | //Event Selection// | |
218 | /////////////////// | |
219 | Double_t deltaTdc = 0; | |
220 | Double_t sumTdc = 0; | |
65abd48b | 221 | |
241cca4d | 222 | for(Int_t i = 0; i < 4; ++i) { |
223 | if(esdZDC->GetZDCTDCData(10,i) != 0) { | |
224 | Double_t tdcCnoCorr = 0.025*(esdZDC->GetZDCTDCData(10,i)- | |
225 | esdZDC->GetZDCTDCData(14,i)); | |
226 | Double_t tdcC = esdZDC->GetZDCTDCCorrected(10,i); | |
227 | for(Int_t j = 0; j < 4; ++j) { | |
228 | if(esdZDC->GetZDCTDCData(12,j) != 0) { | |
229 | Double_t tdcAnoCorr = 0.025*(esdZDC->GetZDCTDCData(12,j)- | |
230 | esdZDC->GetZDCTDCData(14,j)); | |
231 | Double_t tdcA = esdZDC->GetZDCTDCCorrected(12,j); | |
232 | if(esdZDC->TestBit(AliESDZDC::kCorrectedTDCFilled)) { | |
233 | deltaTdc = tdcC-tdcA; | |
234 | sumTdc = tdcC+tdcA; | |
235 | } | |
236 | else { | |
237 | deltaTdc = tdcCnoCorr-tdcAnoCorr; | |
238 | sumTdc = tdcCnoCorr+tdcAnoCorr; | |
239 | } | |
240 | } | |
241 | } | |
242 | } | |
243 | } | |
65abd48b | 244 | |
241cca4d | 245 | /////////////////////// |
246 | //Global Event Filter// | |
247 | /////////////////////// | |
248 | Bool_t zdcAccSat[21]; | |
65abd48b | 249 | |
241cca4d | 250 | // Bool_t zdcAccSatRunClass[21]; |
251 | for(Int_t k = -10; k <= 10; k++) zdcAccSat[k+10] = kFALSE; | |
65abd48b | 252 | |
241cca4d | 253 | |
254 | if(deltaTdc!=0. || sumTdc!=0.) { | |
255 | for (Int_t k = -10; k <= 10; ++k) { | |
256 | Float_t zsat = 2.55F * k; | |
65abd48b | 257 | |
241cca4d | 258 | if(k==0) { |
259 | if(((deltaTdc-kZDCrefDelta)*(deltaTdc-kZDCrefDelta)/ | |
260 | (kZDCsigmaDelta*kZDCsigmaDelta)+ | |
261 | (sumTdc -kZDCrefSum ) * (sumTdc -kZDCrefSum ) / | |
262 | (kZDCsigmaSum *kZDCsigmaSum))<=1.0) { | |
263 | zdcAccSat[k+10] = kTRUE; | |
264 | } | |
265 | } | |
266 | else { | |
267 | if(((deltaTdc-kZDCrefDelta-zsat)*(deltaTdc-kZDCrefDelta-zsat)/ | |
268 | (kZDCsigmaDeltaSat*kZDCsigmaDeltaSat)+ | |
269 | (sumTdc -kZDCrefSum -zsat)*(sumTdc -kZDCrefSum -zsat)/ | |
270 | (kZDCsigmaSumSat *kZDCsigmaSumSat ))<=1.0) { | |
271 | zdcAccSat[k+10] = kTRUE; | |
272 | } | |
273 | } | |
274 | } | |
65abd48b | 275 | } |
276 | ||
241cca4d | 277 | for(Int_t k=-10;k<=10;++k) { |
278 | if(zdcAccSat[k+10] && k!=0 ) { | |
279 | //std::cout<<"Displaced vertex at "<<37.5*(Float_t)k<<" cm"<<std::endl; | |
280 | zvtx = 37.5*(Float_t)k; | |
281 | } | |
282 | } | |
65abd48b | 283 | |
241cca4d | 284 | return zvtx; |
65abd48b | 285 | |
286 | } | |
287 | //____________________________________________________________________ | |
288 | Double_t | |
289 | AliDisplacedVertexSelection::CalculateDisplacedVertexCent(const AliESDEvent* esd) const | |
290 | { | |
291 | Float_t kZDCrefSum =-66.5; | |
292 | Float_t kZDCrefDelta =-2.10; | |
293 | Float_t kZDCsigmaSum = 3.25; | |
294 | Float_t kZDCsigmaDelta = 2.25; | |
295 | Float_t kZDCsigmaSumSat = 2.50; | |
296 | Float_t kZDCsigmaDeltaSat = 1.35; | |
297 | ||
298 | AliESDZDC* esdZDC = esd -> GetESDZDC(); | |
299 | ||
300 | Int_t runnumber = esd->GetRunNumber(); | |
301 | Double_t fcurrentL3 = esd->GetCurrentL3(); | |
302 | Double_t fcurrentDipo = esd->GetCurrentDip(); | |
303 | Double_t cent = -1; | |
304 | ||
241cca4d | 305 | Double_t fzdcEn = (esdZDC->GetZDCN1Energy()+ |
306 | esdZDC->GetZDCP1Energy()+ | |
307 | esdZDC->GetZDCN2Energy()+ | |
308 | esdZDC->GetZDCP2Energy()); | |
309 | Double_t fzemEn = (esdZDC->GetZDCEMEnergy(0)+ | |
310 | esdZDC->GetZDCEMEnergy(1))/8.; | |
65abd48b | 311 | |
312 | ||
241cca4d | 313 | /////////////////// |
314 | //Event Selection// | |
315 | /////////////////// | |
316 | Double_t deltaTdc = 0; | |
317 | Double_t sumTdc = 0; | |
65abd48b | 318 | |
241cca4d | 319 | for(Int_t i = 0; i < 4; ++i) { |
320 | if(esdZDC->GetZDCTDCData(10,i) != 0) { | |
321 | Double_t tdcCnoCorr = 0.025*(esdZDC->GetZDCTDCData(10,i)- | |
322 | esdZDC->GetZDCTDCData(14,i)); | |
323 | Double_t tdcC = esdZDC->GetZDCTDCCorrected(10,i); | |
324 | for(Int_t j = 0; j < 4; ++j) { | |
325 | if(esdZDC->GetZDCTDCData(12,j) != 0) { | |
326 | Double_t tdcAnoCorr = 0.025*(esdZDC->GetZDCTDCData(12,j)- | |
327 | esdZDC->GetZDCTDCData(14,j)); | |
328 | Double_t tdcA = esdZDC->GetZDCTDCCorrected(12,j); | |
329 | if(esdZDC->TestBit(AliESDZDC::kCorrectedTDCFilled)) { | |
330 | deltaTdc = tdcC-tdcA; | |
331 | sumTdc = tdcC+tdcA; | |
332 | } | |
333 | else { | |
334 | deltaTdc = tdcCnoCorr-tdcAnoCorr; | |
335 | sumTdc = tdcCnoCorr+tdcAnoCorr; | |
336 | } | |
337 | } | |
338 | } | |
339 | } | |
340 | } | |
65abd48b | 341 | |
241cca4d | 342 | /////////////////////// |
343 | //Global Event Filter// | |
344 | /////////////////////// | |
345 | Bool_t zdcAccSat[21]; | |
346 | // Bool_t zdcAccSatRunClass[21]; | |
347 | for(Int_t k = -10; k <= 10; k++) zdcAccSat[k+10] = kFALSE; | |
65abd48b | 348 | |
241cca4d | 349 | if(deltaTdc!=0. || sumTdc!=0.) { |
350 | for(Int_t k = -10; k <= 10; ++k) { | |
351 | Float_t zsat = 2.55*(Float_t)k; | |
65abd48b | 352 | |
241cca4d | 353 | if(k==0) { |
354 | if(((deltaTdc-kZDCrefDelta)*(deltaTdc-kZDCrefDelta)/ | |
355 | (kZDCsigmaDelta*kZDCsigmaDelta)+ | |
356 | (sumTdc -kZDCrefSum )*(sumTdc -kZDCrefSum )/ | |
357 | (kZDCsigmaSum *kZDCsigmaSum))<=1.0) { | |
358 | zdcAccSat[k+10] = kTRUE; | |
359 | } | |
360 | } | |
361 | else { | |
362 | if(((deltaTdc-kZDCrefDelta-zsat)*(deltaTdc-kZDCrefDelta-zsat)/ | |
363 | (kZDCsigmaDeltaSat*kZDCsigmaDeltaSat)+ | |
364 | (sumTdc -kZDCrefSum -zsat)*(sumTdc -kZDCrefSum -zsat)/ | |
365 | (kZDCsigmaSumSat *kZDCsigmaSumSat ))<=1.0) { | |
366 | zdcAccSat[k+10] = kTRUE; | |
367 | } | |
368 | ||
369 | } | |
370 | ||
371 | } | |
65abd48b | 372 | } |
373 | ||
241cca4d | 374 | Float_t kZEMcorrPlusPlus[21] = {0.6840,0.7879,0.8722,0.9370,0.9837,1.0137, |
375 | 1.0292,1.0327,1.0271,1.0152,1.0000,0.9844, | |
376 | 0.9714,0.9634,0.9626,0.9708,0.9891,1.0181, | |
377 | 1.0574,1.1060,1.1617}; | |
378 | Float_t kZEMcorrMoinsMoins[21] = {0.7082,0.8012,0.8809,0.9447,0.9916,1.0220, | |
379 | 1.0372,1.0395,1.0318,1.0174,1.0000,0.9830, | |
380 | 0.9699,0.9635,0.9662,0.9794,1.0034,1.0371, | |
381 | 1.0782,1.1224,1.1634}; | |
382 | /////////////////// | |
383 | //ZemEn correction// | |
384 | //////////////////// | |
65abd48b | 385 | |
241cca4d | 386 | for(Int_t k = -10; k <= 10; ++k) { |
387 | if(zdcAccSat[k+10]) { | |
388 | //std::cout<<"Vertex at "<<37.5*(Float_t)k<<"cm from cent"<<std::endl; | |
389 | if(fcurrentDipo>0 && fcurrentL3>0) { | |
390 | fzemEn /= kZEMcorrPlusPlus[k+10]; | |
391 | } | |
392 | if(fcurrentDipo<0 && fcurrentL3<0) { | |
393 | fzemEn /= kZEMcorrMoinsMoins[k+10]; | |
394 | } | |
395 | } | |
396 | } | |
65abd48b | 397 | |
241cca4d | 398 | //////////////////////// |
399 | //Centrality selection// | |
400 | //////////////////////// | |
401 | Double_t fzdcPercentile = -1; | |
402 | Float_t slope; | |
65abd48b | 403 | |
404 | ||
241cca4d | 405 | if(runnumber < 137722 && runnumber >= 137161 ) { |
406 | if(fzemEn > 345.) { | |
407 | slope = (fzdcEn + 16992.)/(fzemEn - 345.); | |
408 | slope += 2.23902e+02; | |
409 | fzdcPercentile = (TMath::ATan(slope) - 1.56667)/9.49434e-05; | |
410 | if (fzdcPercentile<0) fzdcPercentile = 0; | |
411 | } | |
412 | else fzdcPercentile = 100; } | |
413 | else if(runnumber < 139172 && runnumber >= 137722) { | |
414 | if(fzemEn > 295.) { | |
415 | slope = (fzdcEn + 15000.)/(fzemEn - 295.); | |
416 | slope += 2.23660e+02; | |
417 | fzdcPercentile = (TMath::ATan(slope) - 1.56664)/8.99571e-05; | |
418 | if (fzdcPercentile<0) fzdcPercentile = 0; | |
419 | } | |
420 | else fzdcPercentile = 100; | |
421 | } | |
422 | else if(runnumber >= 139172) { | |
423 | if(fzemEn > 345.) { | |
424 | slope = (fzdcEn + 16992.)/(fzemEn - 345.); | |
425 | slope += 2.04789e+02; | |
426 | fzdcPercentile = (TMath::ATan(slope) - 1.56629)/1.02768e-04; | |
427 | if (fzdcPercentile<0) fzdcPercentile = 0; | |
428 | } | |
429 | else fzdcPercentile = 100; | |
430 | } | |
431 | else { | |
432 | if(fzemEn > 345.) { | |
433 | slope = (fzdcEn + 16992.)/(fzemEn - 345.); | |
434 | slope += 2.04789e+02; | |
435 | fzdcPercentile = (TMath::ATan(slope) - 1.56629)/1.02768e-04; | |
436 | if(fzdcPercentile<0) fzdcPercentile = 0; | |
437 | } | |
438 | else fzdcPercentile = 100; | |
439 | } | |
440 | ||
441 | if(fzdcPercentile > 0 && fzdcPercentile <100) { | |
442 | //std::cout<<fzdcPercentile<<std::endl; | |
443 | cent = fzdcPercentile; | |
444 | } | |
445 | return cent; | |
65abd48b | 446 | } |
447 | //____________________________________________________________________ | |
448 | // | |
449 | // EOF | |
450 | // |