coverity + some patch for pp
[u/mrichter/AliRoot.git] / PWG2 / FLOW / AliFlowTools / glauberMC / AliGlauberMC.cxx
CommitLineData
ec852657 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////////////////////////////////////////////////////////////////////////////////
17//
18// Glauber MC implementation
19//
20// origin: PHOBOS experiment
21// alification: Mikolaj Krzewicki, Nikhef, mikolaj.krzewicki@cern.ch
7288f660 22// update: You Zhou, Nikhef, yzhou@nikhef.nl
ec852657 23////////////////////////////////////////////////////////////////////////////////
24
25#include <Riostream.h>
26#include <TMath.h>
27#include <TEllipse.h>
28#include <TRandom.h>
29#include <TNamed.h>
30#include <TObjArray.h>
31#include <TNtuple.h>
32#include <TFile.h>
bb442091 33#include <TTree.h>
ec852657 34
35#include "AliGlauberNucleon.h"
36#include "AliGlauberNucleus.h"
37#include "AliGlauberMC.h"
38
39ClassImp(AliGlauberMC)
40
41//______________________________________________________________________________
42AliGlauberMC::AliGlauberMC(Option_t* NA, Option_t* NB, Double_t xsect) :
566fc3d0 43 TNamed(),
bb442091 44 fANucleus(NA),
45 fBNucleus(NB),
566fc3d0 46 fXSect(xsect),
bb442091 47 fNucleonsA(0),
48 fNucleonsB(0),
566fc3d0 49 fAN(0),
998de9b9 50 fQAN(0),
566fc3d0 51 fBN(0),
998de9b9 52 fQBN(0),
bb442091 53 fnt(0),
54 fMeanX2(0),
55 fMeanY2(0),
56 fMeanXY(0),
2c0ff1e6 57 fMeanX2Parts(0),
58 fMeanY2Parts(0),
59 fMeanXYParts(0),
bb442091 60 fMeanXParts(0),
61 fMeanYParts(0),
2c0ff1e6 62 fMeanOXParts(0),
63 fMeanOYParts(0),
bb442091 64 fMeanXColl(0),
65 fMeanYColl(0),
2c0ff1e6 66 fMeanOXColl(0),
67 fMeanOYColl(0),
bb442091 68 fMeanX2Coll(0),
2c0ff1e6 69 fMeanY2Coll(0),
bb442091 70 fMeanXYColl(0),
4d264dfa 71 fMeanXCom(0),
72 fMeanYCom(0),
73 fMeanOXCom(0),
74 fMeanOYCom(0),
75 fMeanX2Com(0),
76 fMeanY2Com(0),
77 fMeanXYCom(0),
bb442091 78 fMeanXSystem(0),
79 fMeanYSystem(0),
43215add 80 fMeanXA(0),
81 fMeanYA(0),
82 fMeanXB(0),
83 fMeanYB(0),
2c0ff1e6 84 fMeanOXA(0),
85 fMeanOYA(0),
86 fMeanOXB(0),
87 fMeanOYB(0),
43215add 88 fBMC(0),
bb442091 89 fEvents(0),
90 fTotalEvents(0),
566fc3d0 91 fBMin(0.),
92 fBMax(20.),
7288f660 93 fMultType(kNBDSV),
bb442091 94 fMaxNpartFound(0),
2c0ff1e6 95 fONpart(0),
96 fONcoll(0),
4d264dfa 97 fONcom(0),
bb442091 98 fNpart(0),
99 fNcoll(0),
4d264dfa 100 fNcom(0),
106d1ca1 101 fMeanr2(0),
102 fMeanr3(0),
103 fMeanr4(0),
104 fMeanr5(0),
105 fMeanr2Cos2Phi(0),
106 fMeanr2Sin2Phi(0),
107 fMeanr2Cos3Phi(0),
108 fMeanr2Sin3Phi(0),
109 fMeanr2Cos4Phi(0),
110 fMeanr2Sin4Phi(0),
111 fMeanr2Cos5Phi(0),
112 fMeanr2Sin5Phi(0),
113 fMeanr3Cos3Phi(0),
114 fMeanr3Sin3Phi(0),
115 fMeanr4Cos4Phi(0),
116 fMeanr4Sin4Phi(0),
117 fMeanr5Cos5Phi(0),
118 fMeanr5Sin5Phi(0),
9a5780fa 119 fMeanr2Coll(0),
120 fMeanr3Coll(0),
121 fMeanr4Coll(0),
122 fMeanr5Coll(0),
123 fMeanr2Cos2PhiColl(0),
124 fMeanr2Sin2PhiColl(0),
125 fMeanr2Cos3PhiColl(0),
126 fMeanr2Sin3PhiColl(0),
127 fMeanr2Cos4PhiColl(0),
128 fMeanr2Sin4PhiColl(0),
129 fMeanr2Cos5PhiColl(0),
130 fMeanr2Sin5PhiColl(0),
131 fMeanr3Cos3PhiColl(0),
132 fMeanr3Sin3PhiColl(0),
133 fMeanr4Cos4PhiColl(0),
134 fMeanr4Sin4PhiColl(0),
135 fMeanr5Cos5PhiColl(0),
136 fMeanr5Sin5PhiColl(0),
4d264dfa 137 fMeanr2Com(0),
138 fMeanr3Com(0),
139 fMeanr4Com(0),
140 fMeanr5Com(0),
141 fMeanr2Cos2PhiCom(0),
142 fMeanr2Sin2PhiCom(0),
143 fMeanr2Cos3PhiCom(0),
144 fMeanr2Sin3PhiCom(0),
145 fMeanr2Cos4PhiCom(0),
146 fMeanr2Sin4PhiCom(0),
147 fMeanr2Cos5PhiCom(0),
148 fMeanr2Sin5PhiCom(0),
149 fMeanr3Cos3PhiCom(0),
150 fMeanr3Sin3PhiCom(0),
151 fMeanr4Cos4PhiCom(0),
152 fMeanr4Sin4PhiCom(0),
153 fMeanr5Cos5PhiCom(0),
154 fMeanr5Sin5PhiCom(0),
2c0ff1e6 155 fSx2Parts(0.),
156 fSy2Parts(0.),
157 fSxyParts(0.),
566fc3d0 158 fSx2Coll(0.),
159 fSy2Coll(0.),
160 fSxyColl(0.),
4d264dfa 161 fSx2Com(0.),
162 fSy2Com(0.),
163 fSxyCom(0.),
bb442091 164 fX(0.13),
14244e14 165 fNpp(8.),
fee95d0f 166 fDoPartProd(kFALSE)
ec852657 167{
43215add 168 //ctor
7288f660 169 for (UInt_t i=0; i<(sizeof(fdNdEtaParam)/sizeof(fdNdEtaParam[0])); i++)
170 {
171 fdNdEtaParam[i]=0.0;
172 }
bb442091 173
566fc3d0 174 SetName(Form("Glauber_%s_%s",fANucleus.GetName(),fBNucleus.GetName()));
175 SetTitle(Form("Glauber %s+%s Version",fANucleus.GetName(),fBNucleus.GetName()));
176}
177
178//______________________________________________________________________________
179AliGlauberMC::~AliGlauberMC()
180{
181 //dtor
182 delete fnt;
183}
184
185//______________________________________________________________________________
186AliGlauberMC::AliGlauberMC(const AliGlauberMC& in):
187 TNamed(in),
188 fANucleus(in.fANucleus),
189 fBNucleus(in.fBNucleus),
190 fXSect(in.fXSect),
191 fNucleonsA(in.fNucleonsA),
192 fNucleonsB(in.fNucleonsB),
193 fAN(in.fAN),
998de9b9 194 fQAN(in.fQAN),
566fc3d0 195 fBN(in.fBN),
998de9b9 196 fQBN(in.fQBN),
566fc3d0 197 fnt(in.fnt),
198 fMeanX2(in.fMeanX2),
199 fMeanY2(in.fMeanY2),
200 fMeanXY(in.fMeanXY),
2c0ff1e6 201 fMeanX2Parts(in.fMeanX2Parts),
202 fMeanY2Parts(in.fMeanY2Parts),
203 fMeanXYParts(in.fMeanXYParts),
566fc3d0 204 fMeanXParts(in.fMeanXParts),
205 fMeanYParts(in.fMeanYParts),
2c0ff1e6 206 fMeanOXParts(in.fMeanOXParts),
207 fMeanOYParts(in.fMeanOYParts),
566fc3d0 208 fMeanXColl(in.fMeanXColl),
209 fMeanYColl(in.fMeanYColl),
2c0ff1e6 210 fMeanOXColl(in.fMeanOXColl),
211 fMeanOYColl(in.fMeanOYColl),
566fc3d0 212 fMeanX2Coll(in.fMeanX2Coll),
213 fMeanY2Coll(in.fMeanY2Coll),
214 fMeanXYColl(in.fMeanXYColl),
4d264dfa 215 fMeanXCom(in.fMeanXCom),
216 fMeanYCom(in.fMeanYCom),
217 fMeanOXCom(in.fMeanOXCom),
218 fMeanOYCom(in.fMeanOYCom),
219 fMeanX2Com(in.fMeanX2Com),
220 fMeanY2Com(in.fMeanY2Com),
221 fMeanXYCom(in.fMeanXYCom),
566fc3d0 222 fMeanXSystem(in.fMeanXSystem),
223 fMeanYSystem(in.fMeanYSystem),
43215add 224 fMeanXA(in.fMeanXA),
225 fMeanYA(in.fMeanYA),
226 fMeanXB(in.fMeanXB),
227 fMeanYB(in.fMeanYB),
2c0ff1e6 228 fMeanOXA(in.fMeanOXA),
229 fMeanOYA(in.fMeanOYA),
230 fMeanOXB(in.fMeanOXB),
231 fMeanOYB(in.fMeanOYB),
43215add 232 fBMC(in.fBMC),
566fc3d0 233 fEvents(in.fEvents),
234 fTotalEvents(in.fTotalEvents),
235 fBMin(in.fBMin),
236 fBMax(in.fBMax),
7288f660 237 fMultType(in.fMultType),
566fc3d0 238 fMaxNpartFound(in.fMaxNpartFound),
2c0ff1e6 239 fONpart(in.fONpart),
240 fONcoll(in.fONcoll),
4d264dfa 241 fONcom(in.fONcom),
566fc3d0 242 fNpart(in.fNpart),
243 fNcoll(in.fNcoll),
4d264dfa 244 fNcom(in.fNcom),
106d1ca1 245 fMeanr2(in.fMeanr2),
246 fMeanr3(in.fMeanr3),
247 fMeanr4(in.fMeanr4),
248 fMeanr5(in.fMeanr5),
249 fMeanr2Cos2Phi(in.fMeanr2Cos2Phi),
250 fMeanr2Sin2Phi(in.fMeanr2Sin2Phi),
251 fMeanr2Cos3Phi(in.fMeanr2Cos3Phi),
252 fMeanr2Sin3Phi(in.fMeanr2Sin3Phi),
253 fMeanr2Cos4Phi(in.fMeanr2Cos4Phi),
254 fMeanr2Sin4Phi(in.fMeanr2Sin4Phi),
255 fMeanr2Cos5Phi(in.fMeanr2Cos5Phi),
256 fMeanr2Sin5Phi(in.fMeanr2Sin5Phi),
257 fMeanr3Cos3Phi(in.fMeanr3Cos3Phi),
258 fMeanr3Sin3Phi(in.fMeanr3Sin3Phi),
259 fMeanr4Cos4Phi(in.fMeanr4Cos4Phi),
260 fMeanr4Sin4Phi(in.fMeanr4Sin4Phi),
261 fMeanr5Cos5Phi(in.fMeanr5Cos5Phi),
262 fMeanr5Sin5Phi(in.fMeanr5Sin5Phi),
9a5780fa 263 fMeanr2Coll(in.fMeanr2Coll),
264 fMeanr3Coll(in.fMeanr3Coll),
265 fMeanr4Coll(in.fMeanr4Coll),
266 fMeanr5Coll(in.fMeanr5Coll),
267 fMeanr2Cos2PhiColl(in.fMeanr2Cos2PhiColl),
268 fMeanr2Sin2PhiColl(in.fMeanr2Sin2PhiColl),
269 fMeanr2Cos3PhiColl(in.fMeanr2Cos3PhiColl),
270 fMeanr2Sin3PhiColl(in.fMeanr2Sin3PhiColl),
271 fMeanr2Cos4PhiColl(in.fMeanr2Cos4PhiColl),
272 fMeanr2Sin4PhiColl(in.fMeanr2Sin4PhiColl),
273 fMeanr2Cos5PhiColl(in.fMeanr2Cos5PhiColl),
274 fMeanr2Sin5PhiColl(in.fMeanr2Sin5PhiColl),
275 fMeanr3Cos3PhiColl(in.fMeanr3Cos3PhiColl),
276 fMeanr3Sin3PhiColl(in.fMeanr3Sin3PhiColl),
277 fMeanr4Cos4PhiColl(in.fMeanr4Cos4PhiColl),
278 fMeanr4Sin4PhiColl(in.fMeanr4Sin4PhiColl),
279 fMeanr5Cos5PhiColl(in.fMeanr5Cos5PhiColl),
280 fMeanr5Sin5PhiColl(in.fMeanr5Sin5PhiColl),
4d264dfa 281 fMeanr2Com(in.fMeanr2Com),
282 fMeanr3Com(in.fMeanr3Com),
283 fMeanr4Com(in.fMeanr4Com),
284 fMeanr5Com(in.fMeanr5Com),
285 fMeanr2Cos2PhiCom(in.fMeanr2Cos2PhiCom),
286 fMeanr2Sin2PhiCom(in.fMeanr2Sin2PhiCom),
287 fMeanr2Cos3PhiCom(in.fMeanr2Cos3PhiCom),
288 fMeanr2Sin3PhiCom(in.fMeanr2Sin3PhiCom),
289 fMeanr2Cos4PhiCom(in.fMeanr2Cos4PhiCom),
290 fMeanr2Sin4PhiCom(in.fMeanr2Sin4PhiCom),
291 fMeanr2Cos5PhiCom(in.fMeanr2Cos5PhiCom),
292 fMeanr2Sin5PhiCom(in.fMeanr2Sin5PhiCom),
293 fMeanr3Cos3PhiCom(in.fMeanr3Cos3PhiCom),
294 fMeanr3Sin3PhiCom(in.fMeanr3Sin3PhiCom),
295 fMeanr4Cos4PhiCom(in.fMeanr4Cos4PhiCom),
296 fMeanr4Sin4PhiCom(in.fMeanr4Sin4PhiCom),
297 fMeanr5Cos5PhiCom(in.fMeanr5Cos5PhiCom),
298 fMeanr5Sin5PhiCom(in.fMeanr5Sin5PhiCom),
2c0ff1e6 299 fSx2Parts(in.fSx2Parts),
300 fSy2Parts(in.fSy2Parts),
301 fSxyParts(in.fSxyParts),
566fc3d0 302 fSx2Coll(in.fSx2Coll),
303 fSy2Coll(in.fSy2Coll),
304 fSxyColl(in.fSxyColl),
4d264dfa 305 fSx2Com(in.fSx2Com),
306 fSy2Com(in.fSy2Com),
307 fSxyCom(in.fSxyCom),
566fc3d0 308 fX(in.fX),
97d98e48 309 fNpp(in.fNpp),
fee95d0f 310 fDoPartProd(kFALSE)
566fc3d0 311{
312 //copy ctor
7288f660 313 memcpy(fdNdEtaParam,in.fdNdEtaParam,sizeof(fdNdEtaParam));
566fc3d0 314}
315
316//______________________________________________________________________________
317AliGlauberMC& AliGlauberMC::operator=(const AliGlauberMC& in)
318{
319 //assignment
2047680a 320 if (&in==this) return *this;
7288f660 321 fANucleus=in.fANucleus;
322 fBNucleus=in.fBNucleus;
323 fXSect=in.fXSect;
566fc3d0 324 fNucleonsA=in.fNucleonsA;
325 fNucleonsB=in.fNucleonsB;
7288f660 326 fAN=in.fAN;
998de9b9 327 fQAN=in.fQAN;
7288f660 328 fBN=in.fBN;
998de9b9 329 fQBN=in.fQBN;
7288f660 330 fnt=in.fnt;
331 fMeanX2=in.fMeanX2;
332 fMeanY2=in.fMeanY2;
333 fMeanXY=in.fMeanXY;
106d1ca1 334 fMeanr2=in.fMeanr2;
335 fMeanr3=in.fMeanr3;
336 fMeanr4=in.fMeanr4;
337 fMeanr5=in.fMeanr5;
566fc3d0 338 fMeanXParts=in.fMeanXParts;
339 fMeanYParts=in.fMeanYParts;
2c0ff1e6 340 fMeanOXParts=in.fMeanOXParts;
341 fMeanOYParts=in.fMeanOYParts;
566fc3d0 342 fMeanXColl=in.fMeanXColl;
343 fMeanYColl=in.fMeanYColl;
2c0ff1e6 344 fMeanOXColl=in.fMeanOXColl;
345 fMeanOYColl=in.fMeanOYColl;
566fc3d0 346 fMeanX2Coll=in.fMeanX2Coll;
347 fMeanY2Coll=in.fMeanY2Coll;
348 fMeanXYColl=in.fMeanXYColl;
9a5780fa 349 fMeanr2Coll=in.fMeanr2Coll;
350 fMeanr3Coll=in.fMeanr3Coll;
351 fMeanr4Coll=in.fMeanr4Coll;
352 fMeanr5Coll=in.fMeanr5Coll;
4d264dfa 353 fMeanXCom=in.fMeanXColl;
354 fMeanYCom=in.fMeanYColl;
355 fMeanOXCom=in.fMeanOXCom;
356 fMeanOYCom=in.fMeanOYCom;
357 fMeanX2Com=in.fMeanX2Com;
358 fMeanY2Com=in.fMeanY2Com;
359 fMeanXYCom=in.fMeanXYCom;
360 fMeanr2Com=in.fMeanr2Com;
361 fMeanr3Com=in.fMeanr3Com;
362 fMeanr4Com=in.fMeanr4Com;
363 fMeanr5Com=in.fMeanr5Com;
566fc3d0 364 fMeanXSystem=in.fMeanXSystem;
2c0ff1e6 365 fMeanYSystem=in.fMeanYSystem;
7288f660 366 fMeanXA=in.fMeanXA;
367 fMeanYA=in.fMeanYA;
368 fMeanXB=in.fMeanXB;
369 fMeanYB=in.fMeanYB;
2c0ff1e6 370 fMeanOXA=in.fMeanOXA;
371 fMeanOYA=in.fMeanOYA;
372 fMeanOXB=in.fMeanOXB;
373 fMeanOYB=in.fMeanOYB;
7288f660 374 fBMC=in.fBMC;
375 fEvents=in.fEvents;
566fc3d0 376 fTotalEvents=in.fTotalEvents;
7288f660 377 fBMin=in.fBMin;
378 fBMax=in.fBMax;
379 fMultType=in.fMultType,
380 memcpy(fdNdEtaParam,in.fdNdEtaParam,sizeof(fdNdEtaParam));
566fc3d0 381 fMaxNpartFound=in.fMaxNpartFound;
7288f660 382 fNpart=in.fNpart;
383 fNcoll=in.fNcoll;
4d264dfa 384 fNcom=in.fNcom;
2c0ff1e6 385 fONpart=in.fONpart;
386 fONcoll=in.fONcoll;
4d264dfa 387 fONcom=in.fONcom;
106d1ca1 388 fMeanr2Cos2Phi=in.fMeanr2Cos2Phi;
389 fMeanr2Sin2Phi=in.fMeanr2Sin2Phi;
390 fMeanr2Cos3Phi=in.fMeanr2Cos3Phi;
391 fMeanr2Sin3Phi=in.fMeanr2Sin3Phi;
392 fMeanr2Cos4Phi=in.fMeanr2Cos4Phi;
393 fMeanr2Sin4Phi=in.fMeanr2Sin4Phi;
394 fMeanr2Cos5Phi=in.fMeanr2Cos5Phi;
395 fMeanr2Sin5Phi=in.fMeanr2Sin5Phi;
396 fMeanr3Cos3Phi=in.fMeanr3Cos3Phi;
397 fMeanr3Sin3Phi=in.fMeanr3Sin3Phi;
398 fMeanr4Cos4Phi=in.fMeanr4Cos4Phi;
399 fMeanr4Sin4Phi=in.fMeanr4Sin4Phi;
400 fMeanr5Cos5Phi=in.fMeanr5Cos5Phi;
401 fMeanr5Sin5Phi=in.fMeanr5Sin5Phi;
9a5780fa 402 fMeanr2Cos2PhiColl=in.fMeanr2Cos2PhiColl;
403 fMeanr2Sin2PhiColl=in.fMeanr2Sin2PhiColl;
404 fMeanr2Cos3PhiColl=in.fMeanr2Cos3PhiColl;
405 fMeanr2Sin3PhiColl=in.fMeanr2Sin3PhiColl;
406 fMeanr2Cos4PhiColl=in.fMeanr2Cos4PhiColl;
407 fMeanr2Sin4PhiColl=in.fMeanr2Sin4PhiColl;
408 fMeanr2Cos5PhiColl=in.fMeanr2Cos5PhiColl;
409 fMeanr2Sin5PhiColl=in.fMeanr2Sin5PhiColl;
410 fMeanr3Cos3PhiColl=in.fMeanr3Cos3PhiColl;
411 fMeanr3Sin3PhiColl=in.fMeanr3Sin3PhiColl;
412 fMeanr4Cos4PhiColl=in.fMeanr4Cos4PhiColl;
413 fMeanr4Sin4PhiColl=in.fMeanr4Sin4PhiColl;
414 fMeanr5Cos5PhiColl=in.fMeanr5Cos5PhiColl;
415 fMeanr5Sin5PhiColl=in.fMeanr5Sin5PhiColl;
4d264dfa 416 fMeanr2Cos2PhiCom=in.fMeanr2Cos2PhiCom;
417 fMeanr2Sin2PhiCom=in.fMeanr2Sin2PhiCom;
418 fMeanr2Cos3PhiCom=in.fMeanr2Cos3PhiCom;
419 fMeanr2Sin3PhiCom=in.fMeanr2Sin3PhiCom;
420 fMeanr2Cos4PhiCom=in.fMeanr2Cos4PhiCom;
421 fMeanr2Sin4PhiCom=in.fMeanr2Sin4PhiCom;
422 fMeanr2Cos5PhiCom=in.fMeanr2Cos5PhiCom;
423 fMeanr2Sin5PhiCom=in.fMeanr2Sin5PhiCom;
424 fMeanr3Cos3PhiCom=in.fMeanr3Cos3PhiCom;
425 fMeanr3Sin3PhiCom=in.fMeanr3Sin3PhiCom;
426 fMeanr4Cos4PhiCom=in.fMeanr4Cos4PhiCom;
427 fMeanr4Sin4PhiCom=in.fMeanr4Sin4PhiCom;
428 fMeanr5Cos5PhiCom=in.fMeanr5Cos5PhiCom;
429 fMeanr5Sin5PhiCom=in.fMeanr5Sin5PhiCom;
2c0ff1e6 430 fSx2Parts=in.fSx2Parts;
431 fSy2Parts=in.fSy2Parts;
432 fSxyParts=in.fSxyParts;
7288f660 433 fSx2Coll=in.fSx2Coll;
434 fSy2Coll=in.fSy2Coll;
435 fSxyColl=in.fSxyColl;
4d264dfa 436 fSx2Com=in.fSx2Com;
437 fSy2Com=in.fSy2Com;
438 fSxyCom=in.fSxyCom;
7288f660 439 fX=in.fX;
440 fNpp=in.fNpp;
566fc3d0 441 return *this;
ec852657 442}
443
444//______________________________________________________________________________
445Bool_t AliGlauberMC::CalcEvent(Double_t bgen)
446{
bb442091 447 // prepare event
448 fANucleus.ThrowNucleons(-bgen/2.);
449 fNucleonsA = fANucleus.GetNucleons();
450 fAN = fANucleus.GetN();
998de9b9 451 fQAN = fAN * 3;
452 //fAN = 3 * fANucleus.GetN(); // for Pb, Number of quark = 3*208;
bb442091 453 for (Int_t i = 0; i<fAN; i++)
7288f660 454 {
455 AliGlauberNucleon *nucleonA=(AliGlauberNucleon*)(fNucleonsA->UncheckedAt(i));
456 nucleonA->SetInNucleusA();
457 }
bb442091 458 fBNucleus.ThrowNucleons(bgen/2.);
459 fNucleonsB = fBNucleus.GetNucleons();
998de9b9 460 //fBN = 3 * fBNucleus.GetN(); // Number of quark = number of nucleus*3;
bb442091 461 fBN = fBNucleus.GetN();
998de9b9 462 fQBN = fBN * 3;
bb442091 463 for (Int_t i = 0; i<fBN; i++)
7288f660 464 {
465 AliGlauberNucleon *nucleonB=(AliGlauberNucleon*)(fNucleonsB->UncheckedAt(i));
466 nucleonB->SetInNucleusB();
467 }
468
bb442091 469 // "ball" diameter = distance at which two balls interact
470 Double_t d2 = (Double_t)fXSect/(TMath::Pi()*10); // in fm^2
7288f660 471
bb442091 472 // for each of the A nucleons in nucleus B
473 for (Int_t i = 0; i<fBN; i++)
7288f660 474 {
475 AliGlauberNucleon *nucleonB=(AliGlauberNucleon*)(fNucleonsB->UncheckedAt(i));
476 for (Int_t j = 0 ; j < fAN ; j++)
bb442091 477 {
7288f660 478 AliGlauberNucleon *nucleonA=(AliGlauberNucleon*)(fNucleonsA->UncheckedAt(j));
479 Double_t dx = nucleonB->GetX()-nucleonA->GetX();
480 Double_t dy = nucleonB->GetY()-nucleonA->GetY();
481 Double_t dij = dx*dx+dy*dy;
482 if (dij < d2)
483 {
484 nucleonB->Collide();
485 nucleonA->Collide();
486 }
bb442091 487 }
7288f660 488 }
489
ec852657 490 return CalcResults(bgen);
491}
492
493//______________________________________________________________________________
494Bool_t AliGlauberMC::CalcResults(Double_t bgen)
495{
bb442091 496 // calc results for the given event
497 //return true if we have participants
7288f660 498
f98ea3f1 499 fNpart=0;
500 fNcoll=0;
4d264dfa 501 fNcom=0;
2c0ff1e6 502 fONpart=0;
503 fONcoll=0;
4d264dfa 504 fONcom=0;
106d1ca1 505 fMeanX2=0.;
506 fMeanY2=0.;
507 fMeanXY=0.;
508 fMeanXParts=0.;
509 fMeanYParts=0.;
2c0ff1e6 510 fMeanOXParts=0.;
511 fMeanOYParts=0.;
106d1ca1 512 fMeanXColl=0.;
513 fMeanYColl=0.;
2c0ff1e6 514 fMeanOXColl=0.;
515 fMeanOYColl=0.;
106d1ca1 516 fMeanX2Coll=0.;
517 fMeanY2Coll=0.;
518 fMeanXYColl=0.;
4d264dfa 519 fMeanXCom=0.;
520 fMeanYCom=0.;
521 fMeanOXCom=0.;
522 fMeanOYCom=0.;
523 fMeanX2Com=0.;
524 fMeanY2Com=0.;
525 fMeanXYCom=0.;
106d1ca1 526 fMeanXSystem=0.;
527 fMeanYSystem=0.;
528 fMeanXA=0.;
529 fMeanYA=0.;
530 fMeanXB=0.;
531 fMeanYB=0.;
2c0ff1e6 532 fMeanOXA=0.;
533 fMeanOYA=0.;
534 fMeanOXB=0.;
535 fMeanOYB=0.;
106d1ca1 536 fMeanr2=0.;
537 fMeanr3=0.;
538 fMeanr4=0.;
539 fMeanr5=0.;
540 fMeanr2Cos2Phi=0.;
541 fMeanr2Sin2Phi=0.;
542 fMeanr2Cos3Phi=0.;
543 fMeanr2Sin3Phi=0.;
544 fMeanr2Cos4Phi=0.;
545 fMeanr2Sin4Phi=0.;
546 fMeanr2Cos5Phi=0.;
547 fMeanr2Sin5Phi=0.;
548 fMeanr3Cos3Phi=0.;
549 fMeanr3Sin3Phi=0.;
550 fMeanr4Cos4Phi=0.;
551 fMeanr4Sin4Phi=0.;
552 fMeanr5Cos5Phi=0.;
553 fMeanr5Sin5Phi=0.;
9a5780fa 554 fMeanr2Coll=0.;
555 fMeanr3Coll=0.;
556 fMeanr4Coll=0.;
557 fMeanr5Coll=0.;
558 fMeanr2Cos2PhiColl=0.;
559 fMeanr2Sin2PhiColl=0.;
560 fMeanr2Cos3PhiColl=0.;
561 fMeanr2Sin3PhiColl=0.;
562 fMeanr2Cos4PhiColl=0.;
563 fMeanr2Sin4PhiColl=0.;
564 fMeanr2Cos5PhiColl=0.;
565 fMeanr2Sin5PhiColl=0.;
566 fMeanr3Cos3PhiColl=0.;
567 fMeanr3Sin3PhiColl=0.;
568 fMeanr4Cos4PhiColl=0.;
569 fMeanr4Sin4PhiColl=0.;
570 fMeanr5Cos5PhiColl=0.;
571 fMeanr5Sin5PhiColl=0.;
4d264dfa 572 fMeanr2Com=0.;
573 fMeanr3Com=0.;
574 fMeanr4Com=0.;
575 fMeanr5Com=0.;
576 fMeanr2Cos2PhiCom=0.;
577 fMeanr2Sin2PhiCom=0.;
578 fMeanr2Cos3PhiCom=0.;
579 fMeanr2Sin3PhiCom=0.;
580 fMeanr2Cos4PhiCom=0.;
581 fMeanr2Sin4PhiCom=0.;
582 fMeanr2Cos5PhiCom=0.;
583 fMeanr2Sin5PhiCom=0.;
584 fMeanr3Cos3PhiCom=0.;
585 fMeanr3Sin3PhiCom=0.;
586 fMeanr4Cos4PhiCom=0.;
587 fMeanr4Sin4PhiCom=0.;
588 fMeanr5Cos5PhiCom=0.;
589 fMeanr5Sin5PhiCom=0.;
7288f660 590
bb442091 591 for (Int_t i = 0; i<fAN; i++)
7288f660 592 {
593 AliGlauberNucleon *nucleonA=(AliGlauberNucleon*)(fNucleonsA->UncheckedAt(i));
2c0ff1e6 594 Double_t oXA = nucleonA->GetX();
595 Double_t oYA = nucleonA->GetY();
596 //fMeanOXSystem += oXA;
597 //fMeanOYSystem += oYA;
598 fMeanOXA += oXA;
599 fMeanOYA += oYA;
600
601 if(nucleonA->IsWounded())
602 {
603 fONpart++;
604 fMeanOXParts += oXA;
605 fMeanOYParts += oYA;
998de9b9 606
607 fONcom += (1-0.150);
608 fMeanOXCom += oXA*(1-0.150);
609 fMeanOYCom += oXA*(1-0.150);
2c0ff1e6 610 }
611 }
612
613 for (Int_t i = 0; i<fBN; i++)
614 {
615 AliGlauberNucleon *nucleonB=(AliGlauberNucleon*)(fNucleonsB->UncheckedAt(i));
616 Double_t oXB=nucleonB->GetX();
617 Double_t oYB=nucleonB->GetY();
618
619 if(nucleonB->IsWounded())
620 {
621 Int_t oNcoll = nucleonB->GetNColl();
622 fONpart++;
623 fMeanOXParts += oXB;
624 fMeanOXColl += oXB*oNcoll;
998de9b9 625 fMeanOXCom += oXB*((1-0.150)+0.150*oNcoll);
2c0ff1e6 626 fMeanOYParts += oYB;
627 fMeanOYColl += oYB*oNcoll;
998de9b9 628 fMeanOYCom += oYB*((1-0.150)+0.150*oNcoll);
2c0ff1e6 629 fONcoll += oNcoll;
998de9b9 630 fONcom += ((1-0.150)+0.150*oNcoll);
2c0ff1e6 631 }
632 }
633
634 if (fONpart>0)
635 {
636 fMeanOXParts /= fONpart;
637 fMeanOYParts /= fONpart;
638 }
639 else
640 {
641 fMeanOXParts = 0;
642 fMeanOYParts = 0;
643 }
644
645 if (fONcoll>0)
646 {
647 fMeanOXColl /= fONcoll;
648 fMeanOYColl /= fONcoll;
649 }
650 else
651 {
652 fMeanOXColl = 0;
653 fMeanOYColl = 0;
654 }
4d264dfa 655
656 if (fONcom>0)
657 {
658 fMeanOXCom /= fONcom;
659 fMeanOYCom /= fONcom;
660 }
661 else
662 {
663 fMeanOXCom = 0;
664 fMeanOYCom = 0;
665 }
2c0ff1e6 666
667 //////////////////////////////////////////////////////////////////
668 for (Int_t i = 0; i<fAN; i++)
669 {
670 AliGlauberNucleon *nucleonA=(AliGlauberNucleon*)(fNucleonsA->UncheckedAt(i));
671 Double_t xAA = nucleonA->GetX(); // X
672 Double_t yAA = nucleonA->GetY(); // Y
ddf04a47 673 Double_t xAPart = xAA - fMeanOXParts; // X'
674 Double_t yAPart = yAA - fMeanOYParts; // Y'
675 Double_t r2APart = xAPart *xAPart+yAPart*yAPart; // r'^2
676 Double_t rAPart = TMath::Sqrt(r2APart); // r'
677 Double_t r3APart = r2APart*rAPart;
678 Double_t r4APart = r3APart*rAPart;
679 Double_t r5APart = r4APart*rAPart;
680 Double_t phiAPart = TMath::ATan2(yAPart,xAPart);
681 Double_t sin2PhiAPart = TMath::Sin(2.*phiAPart);
682 Double_t cos2PhiAPart = TMath::Cos(2.*phiAPart);
683 Double_t sin3PhiAPart = TMath::Sin(3.*phiAPart);
684 Double_t cos3PhiAPart = TMath::Cos(3.*phiAPart);
685 Double_t sin4PhiAPart = TMath::Sin(4.*phiAPart);
686 Double_t cos4PhiAPart = TMath::Cos(4.*phiAPart);
687 Double_t sin5PhiAPart = TMath::Sin(5.*phiAPart);
688 Double_t cos5PhiAPart = TMath::Cos(5.*phiAPart);
689
998de9b9 690 // for combine
691 Double_t xACom = xAA - fMeanOXCom; // X'-Combine
692 Double_t yACom = yAA - fMeanOYCom; // Y'-C
693 Double_t r2ACom = xACom *xACom+yACom*yACom; // r'^2-C
694 Double_t rACom = TMath::Sqrt(r2ACom); // r'-C
695 Double_t r3ACom = r2ACom*rACom;
696 Double_t r4ACom = r3ACom*rACom;
697 Double_t r5ACom = r4ACom*rACom;
698 Double_t phiACom = TMath::ATan2(yACom,xACom);
699 Double_t sin2PhiACom = TMath::Sin(2.*phiACom);
700 Double_t cos2PhiACom = TMath::Cos(2.*phiACom);
701 Double_t sin3PhiACom = TMath::Sin(3.*phiACom);
702 Double_t cos3PhiACom = TMath::Cos(3.*phiACom);
703 Double_t sin4PhiACom = TMath::Sin(4.*phiACom);
704 Double_t cos4PhiACom = TMath::Cos(4.*phiACom);
705 Double_t sin5PhiACom = TMath::Sin(5.*phiACom);
706 Double_t cos5PhiACom = TMath::Cos(5.*phiACom);
707
2c0ff1e6 708 fMeanXSystem += xAA;
709 fMeanYSystem += yAA;
710 fMeanXA += xAA;
711 fMeanYA += yAA;
712 fMeanX2 += xAA * xAA;
713 fMeanY2 += yAA * yAA;
714 fMeanXY += xAA * yAA;
715
7288f660 716 if(nucleonA->IsWounded())
2c0ff1e6 717 {
998de9b9 718 //Wounded
7288f660 719 fNpart++;
ddf04a47 720 fMeanXParts += xAPart;
721 fMeanYParts += yAPart;
722 fMeanX2Parts += xAPart * xAPart;
723 fMeanY2Parts += yAPart * yAPart;
724 fMeanXYParts += xAPart * yAPart;
725 fMeanr2 += r2APart;
726 fMeanr3 += r3APart;
727 fMeanr4 += r4APart;
728 fMeanr5 += r5APart;
729 fMeanr2Cos2Phi += r2APart*cos2PhiAPart;
730 fMeanr2Sin2Phi += r2APart*sin2PhiAPart;
731 fMeanr2Cos3Phi += r2APart*cos3PhiAPart;
732 fMeanr2Sin3Phi += r2APart*sin3PhiAPart;
733 fMeanr2Cos4Phi += r2APart*cos4PhiAPart;
734 fMeanr2Sin4Phi += r2APart*sin4PhiAPart;
735 fMeanr2Cos5Phi += r2APart*cos5PhiAPart;
736 fMeanr2Sin5Phi += r2APart*sin5PhiAPart;
737 fMeanr3Cos3Phi += r3APart*cos3PhiAPart;
738 fMeanr3Sin3Phi += r3APart*sin3PhiAPart;
739 fMeanr4Cos4Phi += r4APart*cos4PhiAPart;
740 fMeanr4Sin4Phi += r4APart*sin4PhiAPart;
741 fMeanr5Cos5Phi += r5APart*cos5PhiAPart;
742 fMeanr5Sin5Phi += r5APart*sin5PhiAPart;
998de9b9 743
744 // Combined
745 fMeanXCom += xACom*(1-0.150);
746 fMeanYCom += yACom*(1-0.150);
747 fMeanX2Com += xACom*xACom*(1-0.150);
748 fMeanY2Com += yACom*yACom*(1-0.150);
749 fMeanXYCom += xACom*yACom*(1-0.150);
750 fNcom += (1-0.150);
751 fMeanr2Com += r2ACom*(1-0.150);
752 fMeanr3Com += r3ACom*(1-0.150);
753 fMeanr4Com += r4ACom*(1-0.150);
754 fMeanr5Com += r5ACom*(1-0.150);
755 fMeanr2Cos2PhiCom += r2ACom*cos2PhiACom*(1-0.150);
756 fMeanr2Sin2PhiCom += r2ACom*sin2PhiACom*(1-0.150);
757 fMeanr2Cos3PhiCom += r2ACom*cos3PhiACom*(1-0.150);
758 fMeanr2Sin3PhiCom += r2ACom*sin3PhiACom*(1-0.150);
759 fMeanr2Cos4PhiCom += r2ACom*cos4PhiACom*(1-0.150);
760 fMeanr2Sin4PhiCom += r2ACom*sin4PhiACom*(1-0.150);
761 fMeanr2Cos5PhiCom += r2ACom*cos5PhiACom*(1-0.150);
762 fMeanr2Sin5PhiCom += r2ACom*sin5PhiACom*(1-0.150);
763 fMeanr3Cos3PhiCom += r3ACom*cos3PhiACom*(1-0.150);
764 fMeanr3Sin3PhiCom += r3ACom*sin3PhiACom*(1-0.150);
765 fMeanr4Cos4PhiCom += r4ACom*cos4PhiACom*(1-0.150);
766 fMeanr4Sin4PhiCom += r4ACom*sin4PhiACom*(1-0.150);
767 fMeanr5Cos5PhiCom += r5ACom*cos5PhiACom*(1-0.150);
768 fMeanr5Sin5PhiCom += r5ACom*sin5PhiACom*(1-0.150);
769
bb442091 770 }
7288f660 771 }
ddf04a47 772
bb442091 773 for (Int_t i = 0; i<fBN; i++)
bb442091 774 {
2c0ff1e6 775 AliGlauberNucleon *nucleonB=(AliGlauberNucleon*)(fNucleonsB->UncheckedAt(i));
776 Double_t xBB = nucleonB->GetX();
777 Double_t yBB = nucleonB->GetY();
ddf04a47 778 // for Wounded
779 Double_t xBPart = xBB - fMeanOXParts; // X'
780 Double_t yBPart = yBB - fMeanOYParts; // Y'
781 Double_t r2BPart = xBPart*xBPart+yBPart*yBPart; // r'^2
782 Double_t rBPart = TMath::Sqrt(r2BPart); // r'
783 Double_t r3BPart = r2BPart*rBPart;
784 Double_t r4BPart = r3BPart*rBPart;
785 Double_t r5BPart = r4BPart*rBPart;
786 Double_t phiBPart = TMath::ATan2(yBPart,xBPart);
787 Double_t sin2PhiBPart = TMath::Sin(2.*phiBPart);
788 Double_t cos2PhiBPart = TMath::Cos(2.*phiBPart);
789 Double_t sin3PhiBPart = TMath::Sin(3.*phiBPart);
790 Double_t cos3PhiBPart = TMath::Cos(3.*phiBPart);
791 Double_t sin4PhiBPart = TMath::Sin(4.*phiBPart);
792 Double_t cos4PhiBPart = TMath::Cos(4.*phiBPart);
793 Double_t sin5PhiBPart = TMath::Sin(5.*phiBPart);
794 Double_t cos5PhiBPart = TMath::Cos(5.*phiBPart);
795 // for Binary
796 Double_t xBColl = xBB - fMeanOXColl; // X'-Binary
797 Double_t yBColl = yBB - fMeanOYColl; // Y'-B
798 Double_t r2BColl = xBColl*xBColl+yBColl*yBColl; // r'^2-B
799 Double_t rBColl = TMath::Sqrt(r2BColl); // r'-B
800 Double_t r3BColl = r2BColl*rBColl;
801 Double_t r4BColl = r3BColl*rBColl;
802 Double_t r5BColl = r4BColl*rBColl;
803 Double_t phiBColl = TMath::ATan2(yBColl,xBColl);
804 Double_t sin2PhiBColl = TMath::Sin(2.*phiBColl);
805 Double_t cos2PhiBColl = TMath::Cos(2.*phiBColl);
806 Double_t sin3PhiBColl = TMath::Sin(3.*phiBColl);
807 Double_t cos3PhiBColl = TMath::Cos(3.*phiBColl);
808 Double_t sin4PhiBColl = TMath::Sin(4.*phiBColl);
809 Double_t cos4PhiBColl = TMath::Cos(4.*phiBColl);
810 Double_t sin5PhiBColl = TMath::Sin(5.*phiBColl);
811 Double_t cos5PhiBColl = TMath::Cos(5.*phiBColl);
812 // for combine
813 Double_t xBCom = xBB - fMeanOXCom; // X'-Combine
814 Double_t yBCom = yBB - fMeanOYCom; // Y'-C
815 Double_t r2BCom = xBCom *xBCom+yBCom*yBCom; // r'^2-C
816 Double_t rBCom = TMath::Sqrt(r2BCom); // r'-C
817 Double_t r3BCom = r2BCom*rBCom;
818 Double_t r4BCom = r3BCom*rBCom;
819 Double_t r5BCom = r4BCom*rBCom;
820 Double_t phiBCom = TMath::ATan2(yBCom,xBCom);
821 Double_t sin2PhiBCom = TMath::Sin(2.*phiBCom);
822 Double_t cos2PhiBCom = TMath::Cos(2.*phiBCom);
823 Double_t sin3PhiBCom = TMath::Sin(3.*phiBCom);
824 Double_t cos3PhiBCom = TMath::Cos(3.*phiBCom);
825 Double_t sin4PhiBCom = TMath::Sin(4.*phiBCom);
826 Double_t cos4PhiBCom = TMath::Cos(4.*phiBCom);
827 Double_t sin5PhiBCom = TMath::Sin(5.*phiBCom);
828 Double_t cos5PhiBCom = TMath::Cos(5.*phiBCom);
829
2c0ff1e6 830 fMeanXSystem += xBB;
831 fMeanYSystem += yBB;
832 fMeanXB += xBB;
833 fMeanYB += yBB;
834 fMeanX2 += xBB*xBB;
835 fMeanY2 += yBB*yBB;
836 fMeanXY += xBB*yBB;
837
838 if(nucleonB->IsWounded())
839 {
840 Int_t ncoll = nucleonB->GetNColl();
841 fNpart++;
ddf04a47 842 fMeanXParts += xBPart;
843 fMeanXColl += xBColl*ncoll;
998de9b9 844 fMeanXCom += xBCom*((1-0.150)+0.150*ncoll);
ddf04a47 845 fMeanYParts += yBPart;
846 fMeanYColl += yBColl*ncoll;
998de9b9 847 fMeanYCom += yBCom*((1-0.150)+0.150*ncoll);
ddf04a47 848 fMeanX2Parts += xBPart * xBPart;
849 fMeanX2Coll += xBColl*xBColl*ncoll;
998de9b9 850 fMeanX2Com += xBCom*xBCom*((1-0.150)+0.150*ncoll);
ddf04a47 851 fMeanY2Parts += yBPart * yBPart;
852 fMeanY2Coll += yBColl*yBColl*ncoll;
998de9b9 853 fMeanY2Com += yBCom*yBCom*((1-0.150)+0.150*ncoll);
ddf04a47 854 fMeanXYParts += xBPart * yBPart;
855 fMeanXYColl += xBColl*yBColl*ncoll;
998de9b9 856 fMeanXYCom += xBCom*yBCom*((1-0.150)+0.150*ncoll);
2c0ff1e6 857 fNcoll += ncoll;
998de9b9 858 fNcom += ((1-0.150)+0.150*ncoll);
ddf04a47 859 fMeanr2 += r2BPart;
860 fMeanr3 += r3BPart;
861 fMeanr4 += r4BPart;
862 fMeanr5 += r5BPart;
863 fMeanr2Cos2Phi += r2BPart*cos2PhiBPart;
864 fMeanr2Sin2Phi += r2BPart*sin2PhiBPart;
865 fMeanr2Cos3Phi += r2BPart*cos3PhiBPart;
866 fMeanr2Sin3Phi += r2BPart*sin3PhiBPart;
867 fMeanr2Cos4Phi += r2BPart*cos4PhiBPart;
868 fMeanr2Sin4Phi += r2BPart*sin4PhiBPart;
869 fMeanr2Cos5Phi += r2BPart*cos5PhiBPart;
870 fMeanr2Sin5Phi += r2BPart*sin5PhiBPart;
871 fMeanr3Cos3Phi += r3BPart*cos3PhiBPart;
872 fMeanr3Sin3Phi += r3BPart*sin3PhiBPart;
873 fMeanr4Cos4Phi += r4BPart*cos4PhiBPart;
874 fMeanr4Sin4Phi += r4BPart*sin4PhiBPart;
875 fMeanr5Cos5Phi += r5BPart*cos5PhiBPart;
876 fMeanr5Sin5Phi += r5BPart*sin5PhiBPart;
877 fMeanr2Coll += r2BColl*ncoll;
878 fMeanr3Coll += r3BColl*ncoll;
879 fMeanr4Coll += r4BColl*ncoll;
880 fMeanr5Coll += r5BColl*ncoll;
881 fMeanr2Cos2PhiColl += r2BColl*cos2PhiBColl*ncoll;
882 fMeanr2Sin2PhiColl += r2BColl*sin2PhiBColl*ncoll;
883 fMeanr2Cos3PhiColl += r2BColl*cos3PhiBColl*ncoll;
884 fMeanr2Sin3PhiColl += r2BColl*sin3PhiBColl*ncoll;
885 fMeanr2Cos4PhiColl += r2BColl*cos4PhiBColl*ncoll;
886 fMeanr2Sin4PhiColl += r2BColl*sin4PhiBColl*ncoll;
887 fMeanr2Cos5PhiColl += r2BColl*cos5PhiBColl*ncoll;
888 fMeanr2Sin5PhiColl += r2BColl*sin5PhiBColl*ncoll;
889 fMeanr3Cos3PhiColl += r3BColl*cos3PhiBColl*ncoll;
890 fMeanr3Sin3PhiColl += r3BColl*sin3PhiBColl*ncoll;
891 fMeanr4Cos4PhiColl += r4BColl*cos4PhiBColl*ncoll;
892 fMeanr4Sin4PhiColl += r4BColl*sin4PhiBColl*ncoll;
893 fMeanr5Cos5PhiColl += r5BColl*cos5PhiBColl*ncoll;
894 fMeanr5Sin5PhiColl += r5BColl*sin5PhiBColl*ncoll;
998de9b9 895 fMeanr2Com += r2BCom*((1-0.150)+0.150*ncoll);
896 fMeanr3Com += r3BCom*((1-0.150)+0.150*ncoll);
897 fMeanr4Com += r4BCom*((1-0.150)+0.150*ncoll);
898 fMeanr5Com += r5BCom*((1-0.150)+0.150*ncoll);
899 fMeanr2Cos2PhiCom += r2BCom*cos2PhiBCom*((1-0.150)+0.150*ncoll);
900 fMeanr2Sin2PhiCom += r2BCom*sin2PhiBCom*((1-0.150)+0.150*ncoll);
901 fMeanr2Cos3PhiCom += r2BCom*cos3PhiBCom*((1-0.150)+0.150*ncoll);
902 fMeanr2Sin3PhiCom += r2BCom*sin3PhiBCom*((1-0.150)+0.150*ncoll);
903 fMeanr2Cos4PhiCom += r2BCom*cos4PhiBCom*((1-0.150)+0.150*ncoll);
904 fMeanr2Sin4PhiCom += r2BCom*sin4PhiBCom*((1-0.150)+0.150*ncoll);
905 fMeanr2Cos5PhiCom += r2BCom*cos5PhiBCom*((1-0.150)+0.150*ncoll);
906 fMeanr2Sin5PhiCom += r2BCom*sin5PhiBCom*((1-0.150)+0.150*ncoll);
907 fMeanr3Cos3PhiCom += r3BCom*cos3PhiBCom*((1-0.150)+0.150*ncoll);
908 fMeanr3Sin3PhiCom += r3BCom*sin3PhiBCom*((1-0.150)+0.150*ncoll);
909 fMeanr4Cos4PhiCom += r4BCom*cos4PhiBCom*((1-0.150)+0.150*ncoll);
910 fMeanr4Sin4PhiCom += r4BCom*sin4PhiBCom*((1-0.150)+0.150*ncoll);
911 fMeanr5Cos5PhiCom += r5BCom*cos5PhiBCom*((1-0.150)+0.150*ncoll);
912 fMeanr5Sin5PhiCom += r5BCom*sin5PhiBCom*((1-0.150)+0.150*ncoll);
2c0ff1e6 913 }
bb442091 914 }
2c0ff1e6 915
bb442091 916 if (fNpart>0)
2c0ff1e6 917 {
918 fMeanXParts /= fNpart;
919 fMeanYParts /= fNpart;
920 fMeanX2Parts /= fNpart;
921 fMeanY2Parts /= fNpart;
922 fMeanXYParts /= fNpart;
923 fMeanr2 /= fNpart;
924 fMeanr3 /= fNpart;
925 fMeanr4 /= fNpart;
926 fMeanr5 /= fNpart;
927 fMeanr2Cos2Phi /= fNpart;
928 fMeanr2Sin2Phi /= fNpart;
929 fMeanr2Cos3Phi /= fNpart;
930 fMeanr2Sin3Phi /= fNpart;
931 fMeanr2Cos4Phi /= fNpart;
932 fMeanr2Sin4Phi /= fNpart;
933 fMeanr2Cos5Phi /= fNpart;
934 fMeanr2Sin5Phi /= fNpart;
935 fMeanr3Cos3Phi /= fNpart;
936 fMeanr3Sin3Phi /= fNpart;
937 fMeanr4Cos4Phi /= fNpart;
938 fMeanr4Sin4Phi /= fNpart;
939 fMeanr5Cos5Phi /= fNpart;
940 fMeanr5Sin5Phi /= fNpart;
941 }
942
bb442091 943 else
2c0ff1e6 944 {
945 fMeanXParts = 0;
946 fMeanYParts = 0;
947 fMeanX2Parts = 0;
948 fMeanY2Parts = 0;
949 fMeanXYParts = 0;
950 fMeanr2 = 0;
951 fMeanr3 = 0;
952 fMeanr4 = 0;
953 fMeanr5 = 0;
954 fMeanr2Cos2Phi = 0;
955 fMeanr2Sin2Phi = 0;
956 fMeanr2Cos3Phi = 0;
957 fMeanr2Sin3Phi = 0;
958 fMeanr2Cos4Phi = 0;
959 fMeanr2Sin4Phi = 0;
960 fMeanr2Cos5Phi = 0;
961 fMeanr2Sin5Phi = 0;
962 fMeanr3Cos3Phi = 0;
963 fMeanr3Sin3Phi = 0;
964 fMeanr4Cos4Phi = 0;
965 fMeanr4Sin4Phi = 0;
966 fMeanr5Cos5Phi = 0;
967 fMeanr5Sin5Phi = 0;
968 }
969
bb442091 970 if (fNcoll>0)
2c0ff1e6 971 {
972 fMeanXColl /= fNcoll;
973 fMeanYColl /= fNcoll;
974 fMeanX2Coll /= fNcoll;
975 fMeanY2Coll /= fNcoll;
976 fMeanXYColl /= fNcoll;
977 fMeanr2Coll /= fNcoll;
978 fMeanr3Coll /= fNcoll;
979 fMeanr4Coll /= fNcoll;
980 fMeanr5Coll /= fNcoll;
981 fMeanr2Cos2PhiColl /= fNcoll;
982 fMeanr2Sin2PhiColl /= fNcoll;
983 fMeanr2Cos3PhiColl /= fNcoll;
984 fMeanr2Sin3PhiColl /= fNcoll;
985 fMeanr2Cos4PhiColl /= fNcoll;
986 fMeanr2Sin4PhiColl /= fNcoll;
987 fMeanr2Cos5PhiColl /= fNcoll;
988 fMeanr2Sin5PhiColl /= fNcoll;
989 fMeanr3Cos3PhiColl /= fNcoll;
990 fMeanr3Sin3PhiColl /= fNcoll;
991 fMeanr4Cos4PhiColl /= fNcoll;
992 fMeanr4Sin4PhiColl /= fNcoll;
993 fMeanr5Cos5PhiColl /= fNcoll;
994 fMeanr5Sin5PhiColl /= fNcoll;
7288f660 995 }
2c0ff1e6 996
bb442091 997 else
7288f660 998 {
999 fMeanXColl = 0;
1000 fMeanYColl = 0;
1001 fMeanX2Coll = 0;
1002 fMeanY2Coll = 0;
1003 fMeanXYColl = 0;
1004 fMeanr2Cos2PhiColl =0;
1005 fMeanr2Sin2PhiColl =0;
1006 fMeanr2Cos3PhiColl =0;
1007 fMeanr2Sin3PhiColl =0;
1008 fMeanr2Cos4PhiColl =0;
1009 fMeanr2Sin4PhiColl =0;
1010 fMeanr2Cos5PhiColl =0;
1011 fMeanr2Sin5PhiColl =0;
1012 fMeanr3Cos3PhiColl =0;
1013 fMeanr3Sin3PhiColl =0;
1014 fMeanr4Cos4PhiColl =0;
1015 fMeanr4Sin4PhiColl =0;
1016 fMeanr5Cos5PhiColl =0;
1017 fMeanr5Sin5PhiColl =0;
1018 }
2c0ff1e6 1019
4d264dfa 1020 if (fNcom>0)
1021 {
1022 fMeanXCom /= fNcom;
1023 fMeanYCom /= fNcom;
1024 fMeanX2Com /= fNcom;
1025 fMeanY2Com /= fNcom;
1026 fMeanXYCom /= fNcom;
1027 fMeanr2Com /= fNcom;
1028 fMeanr3Com /= fNcom;
1029 fMeanr4Com /= fNcom;
1030 fMeanr5Com /= fNcom;
1031 fMeanr2Cos2PhiCom /= fNcom;
1032 fMeanr2Sin2PhiCom /= fNcom;
1033 fMeanr2Cos3PhiCom /= fNcom;
1034 fMeanr2Sin3PhiCom /= fNcom;
1035 fMeanr2Cos4PhiCom /= fNcom;
1036 fMeanr2Sin4PhiCom /= fNcom;
1037 fMeanr2Cos5PhiCom /= fNcom;
1038 fMeanr2Sin5PhiCom /= fNcom;
1039 fMeanr3Cos3PhiCom /= fNcom;
1040 fMeanr3Sin3PhiCom /= fNcom;
1041 fMeanr4Cos4PhiCom /= fNcom;
1042 fMeanr4Sin4PhiCom /= fNcom;
1043 fMeanr5Cos5PhiCom /= fNcom;
1044 fMeanr5Sin5PhiCom /= fNcom;
1045 }
1046
1047 else
1048 {
1049 fMeanXCom = 0;
1050 fMeanYCom = 0;
1051 fMeanX2Com = 0;
1052 fMeanY2Com = 0;
1053 fMeanXYCom = 0;
1054 fMeanr2Cos2PhiCom =0;
1055 fMeanr2Sin2PhiCom =0;
1056 fMeanr2Cos3PhiCom =0;
1057 fMeanr2Sin3PhiCom =0;
1058 fMeanr2Cos4PhiCom =0;
1059 fMeanr2Sin4PhiCom =0;
1060 fMeanr2Cos5PhiCom =0;
1061 fMeanr2Sin5PhiCom =0;
1062 fMeanr3Cos3PhiCom =0;
1063 fMeanr3Sin3PhiCom =0;
1064 fMeanr4Cos4PhiCom =0;
1065 fMeanr4Sin4PhiCom =0;
1066 fMeanr5Cos5PhiCom =0;
1067 fMeanr5Sin5PhiCom =0;
1068 }
1069
bb442091 1070 if(fAN+fBN>0)
7288f660 1071 {
1072 fMeanXSystem /= (fAN + fBN);
1073 fMeanYSystem /= (fAN + fBN);
2c0ff1e6 1074 fMeanX2 /= (fAN + fBN);
1075 fMeanY2 /= (fAN + fBN);
1076 fMeanXY /= (fAN + fBN);
7288f660 1077 }
bb442091 1078 else
7288f660 1079 {
1080 fMeanXSystem = 0;
1081 fMeanYSystem = 0;
2c0ff1e6 1082 fMeanX2 = 0;
1083 fMeanY2 = 0;
1084 fMeanXY = 0;
7288f660 1085 }
bb442091 1086 if(fAN>0)
7288f660 1087 {
1088 fMeanXA /= fAN;
1089 fMeanYA /= fAN;
1090 }
bb442091 1091 else
7288f660 1092 {
1093 fMeanXA = 0;
1094 fMeanYA = 0;
1095 }
1096
bb442091 1097 if(fBN>0)
7288f660 1098 {
1099 fMeanXB /= fBN;
1100 fMeanYB /= fBN;
1101 }
bb442091 1102 else
7288f660 1103 {
1104 fMeanXB = 0;
1105 fMeanYB = 0;
1106 }
2c0ff1e6 1107
1108
7288f660 1109
2c0ff1e6 1110
1111 //////////////////////////////////////////////////////////////////
1112 fSx2Parts=fMeanX2Parts-(fMeanXParts*fMeanXParts);
1113 fSy2Parts=fMeanY2Parts-(fMeanYParts*fMeanYParts);
1114 fSxyParts=fMeanXYParts-fMeanXParts*fMeanYParts;
bb442091 1115 fSx2Coll=fMeanX2Coll-(fMeanXColl*fMeanXColl);
1116 fSy2Coll=fMeanY2Coll-(fMeanYColl*fMeanYColl);
1117 fSxyColl=fMeanXYColl-fMeanXColl*fMeanYColl;
4d264dfa 1118 fSx2Com=fMeanX2Com-(fMeanXCom*fMeanXCom);
1119 fSy2Com=fMeanY2Com-(fMeanYCom*fMeanYCom);
1120 fSxyCom=fMeanXYCom-fMeanXCom*fMeanYCom;
43215add 1121 fBMC = bgen;
bb442091 1122 fTotalEvents++;
1123 if (fNpart>0) fEvents++;
bb442091 1124 if (fNpart==0) return kFALSE;
1125 if (fNpart > fMaxNpartFound) fMaxNpartFound = fNpart;
1126
1127 return kTRUE;
ec852657 1128}
1129
1130//______________________________________________________________________________
1131void AliGlauberMC::Draw(Option_t* /*option*/)
1132{
43215add 1133 //draw method
bb442091 1134 fANucleus.Draw(fXSect, 2);
1135 fBNucleus.Draw(fXSect, 4);
1136
1137 TEllipse e;
1138 e.SetFillColor(0);
1139 e.SetLineColor(1);
1140 e.SetLineStyle(2);
1141 e.SetLineWidth(1);
1142 e.DrawEllipse(GetB()/2,0,fBNucleus.GetR(),fBNucleus.GetR(),0,360,0);
1143 e.DrawEllipse(-GetB()/2,0,fANucleus.GetR(),fANucleus.GetR(),0,360,0);
ec852657 1144}
1145
1146//______________________________________________________________________________
1147Double_t AliGlauberMC::GetTotXSect() const
1148{
43215add 1149 //total xsection
bb442091 1150 return (1.*fEvents/fTotalEvents)*TMath::Pi()*fBMax*fBMax/100;
ec852657 1151}
1152
1153//______________________________________________________________________________
1154Double_t AliGlauberMC::GetTotXSectErr() const
1155{
43215add 1156 //total xsection error
bb442091 1157 return GetTotXSect()/TMath::Sqrt((Double_t)fEvents) *
7288f660 1158 TMath::Sqrt(Double_t(1.-fEvents/fTotalEvents));
ec852657 1159}
1160
1161//______________________________________________________________________________
bb442091 1162TObjArray *AliGlauberMC::GetNucleons()
ec852657 1163{
43215add 1164 //get array of nucleons
bb442091 1165 if(!fNucleonsA || !fNucleonsB) return 0;
1166 fNucleonsA->SetOwner(0);
1167 fNucleonsB->SetOwner(0);
1168 TObjArray *allnucleons=new TObjArray(fAN+fBN);
1169 allnucleons->SetOwner();
1170 for (Int_t i = 0; i<fAN; i++)
7288f660 1171 {
1172 allnucleons->Add(fNucleonsA->UncheckedAt(i));
1173 }
bb442091 1174 for (Int_t i = 0; i<fBN; i++)
7288f660 1175 {
1176 allnucleons->Add(fNucleonsB->UncheckedAt(i));
1177 }
bb442091 1178 return allnucleons;
1179}
9a5780fa 1180
bb442091 1181//______________________________________________________________________________
1182Double_t AliGlauberMC::NegativeBinomialDistribution(Int_t x, Int_t k, Double_t nmean)
1183{
43215add 1184 //produce a number distributed acc. neg.bin.dist
bb442091 1185 if(k<=0)
7288f660 1186 {
1187 cout << "Error, zero or negative K" << endl;
1188 return 0;
1189 }
bb442091 1190 return (TMath::Binomial(x+k-1,x))
7288f660 1191 *TMath::Power(((nmean/Double_t(k))/(1+nmean/Double_t(k))),Double_t(x))
1192 *(1/(TMath::Power((1+nmean/Double_t(k)),Double_t(k))));
bb442091 1193}
9a5780fa 1194
bb442091 1195//______________________________________________________________________________
43215add 1196Int_t AliGlauberMC::NegativeBinomialRandom(Int_t k, Double_t nmean) const
bb442091 1197{
43215add 1198 //return random integer from a Negative Binomial Distribution
bb442091 1199 static const Int_t fMaxPlot = 1000;
1200 Double_t array[fMaxPlot];
1201 array[0] = NegativeBinomialDistribution(0,k,nmean);
1202 for (Int_t i=1; i<fMaxPlot; i++)
7288f660 1203 {
1204 array[i] = NegativeBinomialDistribution(i,k,nmean) + array[i-1];
1205 }
bb442091 1206 Double_t r = gRandom->Uniform(0,1);
1207 return TMath::BinarySearch(fMaxPlot,array,r)+2;
9a5780fa 1208
bb442091 1209}
7288f660 1210
1211//______________________________________________________________________________
1212Int_t AliGlauberMC::NegativeBinomialRandomSV(Double_t k, Double_t nbar) const
1213{
1214 // negative binomial distribution generator, S. Voloshin, 09-May-2007
1215 Double_t sum=0.;
1216 Int_t i=0;
1217 Double_t ran=gRandom->Rndm();
1218 Double_t trm=1./pow(1.+nbar/k,k);
1219 if (trm==0.)
1220 {
1221 cout<<"NBD overflow"<<" nbar="<<nbar<<" k="<<k<<endl;
1222 return -1;
1223 }
1224 for(i=0; i<2000 && sum<ran ; i++)
1225 {
1226 sum += trm;
1227 trm *= (k+i)/(i+1.)*(nbar/(nbar+k));
1228 }
1229 return i-1;
1230}
1231
bb442091 1232//______________________________________________________________________________
1233Int_t AliGlauberMC::DoubleNegativeBinomialRandom( Int_t k,
7288f660 1234 Double_t nmean,
1235 Int_t k2,
1236 Double_t nmean2,
1237 Double_t alpha ) const
bb442091 1238{
1239 //return random integer from a Double Negative Binomial Distribution
1240 static const Int_t fMaxPlot = 1000;
1241 Double_t array[fMaxPlot];
1242 array[0] = alpha*NegativeBinomialDistribution(0,k,nmean)+(1-alpha)*NegativeBinomialDistribution(0,k2,nmean2);
1243 for (Int_t i=1; i<fMaxPlot; i++)
7288f660 1244 {
1245 array[i] = alpha*NegativeBinomialDistribution(i,k,nmean)+(1-alpha)*NegativeBinomialDistribution(i,k2,nmean2) + array[i-1];
1246 }
bb442091 1247 Double_t r = gRandom->Uniform(0,1);
1248 return TMath::BinarySearch(fMaxPlot,array,r)+2;
1249}
1250
bb442091 1251//______________________________________________________________________________
7288f660 1252Double_t AliGlauberMC::GetdNdEta() const
bb442091 1253{
7288f660 1254 switch (fMultType)
1255 {
1256 case kSimple:
1257 return GetdNdEtaSimple(fdNdEtaParam);
1258 case kNBD:
1259 return GetdNdEtaNBD(fdNdEtaParam);
1260 case kNBDSV:
1261 return GetdNdEtaNBDSV(fdNdEtaParam);
1262 case kTwoNBD:
1263 return GetdNdEtaTwoNBD(fdNdEtaParam);
1264 case kGBW:
1265 return GetdNdEtaGBW(fdNdEtaParam);
1266 default:
1267 return 0.0;
1268 }
bb442091 1269}
14244e14 1270
bb442091 1271//______________________________________________________________________________
7288f660 1272Double_t AliGlauberMC::GetdNdEtaSimple(const Double_t* p) const
1b4934d5 1273{
1274 //Get particle density per unit of rapidity
bb442091 1275 //using two component model
1276 //Parameters: npp, x
7288f660 1277 Double_t nnp = p[0]; //=8.0
1278 Double_t x = p[1]; //=0.13
1279
bb442091 1280 return nnp*((1.-x)*fNpart/2.+x*fNcoll);
0448b9f3 1281}
1282
1283//______________________________________________________________________________
7288f660 1284Double_t AliGlauberMC::GetdNdEtaGBW( const Double_t* p ) const
0448b9f3 1285{
1286 //Get particle density per unit of rapidity
1287 //using the GBW model
bb442091 1288 //Parameters: delta, lambda, snn
7288f660 1289 Double_t delta = p[0]; //=0.79
1290 Double_t lambda = p[1]; //=0.288
1291 Double_t snn = p[2]; //=30.25
1292
bb442091 1293 return fNpart*0.47*TMath::Sqrt(TMath::Power(snn,lambda))
7288f660 1294 * TMath::Power(fNpart,(1.-delta)/3./delta);
1b4934d5 1295}
1296
14244e14 1297//_______________________________________________________________________________
7288f660 1298Double_t AliGlauberMC::GetdNdEtaNBDSV ( const Double_t* p ) const
1299{
1300 //Get particle density per unit of rapidity (from Sergei)
1301 Double_t npp = p[0]; //=2.49
1302 Double_t ratioSgm2Mu = p[1]; //=1.7
1303 Double_t xhard = p[2]; //=0.13
1304
1305 double knb=npp/(ratioSgm2Mu-1.); // sgm^2/mu = 1+ mu/k
1306 double scale = (1.-xhard)*fNpart/2.+xhard*fNcoll; // effectively get number of pp collisions
1307 float nch1=-99.;
1308 if (knb*scale <1000.)
1309 {
1310 nch1=(float) NegativeBinomialRandomSV( npp*scale, knb*scale );
1311 }
1312 else
1313 {
1314 nch1=(float) NegativeBinomialRandomSV( npp*scale/2., knb*scale/2. ) +
1315 (float) NegativeBinomialRandomSV( npp*scale/2., knb*scale/2. );
1316 }
1317 return nch1;
1318}
1319
1320//_______________________________________________________________________________
1321Double_t AliGlauberMC::GetdNdEtaNBD ( const Double_t* p ) const
bb442091 1322{
1323 //Get particle density per unit of rapidity
1324 //using a aandomized number from a negative binomial distrubution
7288f660 1325 //Parameters: k = related to distribition width=3
1326 // nmean = mean of distribution=4
1327 // beta = set contribution of participants / binary collisions to multiplicity=0.13
1328 Int_t k = TMath::Nint(p[0]);
1329 Double_t nmean = p[1];
1330 Double_t beta = p[2];
1331
bb442091 1332 Double_t mulnp=0.0;
1333 for(Int_t i = 0; i<fNpart; i++)
7288f660 1334 {
1335 mulnp+=NegativeBinomialRandom(k,nmean);
1336 }
bb442091 1337 Double_t mulnb=0.0;
1338 for(Int_t i = 0; i<fNcoll; i++)
7288f660 1339 {
1340 mulnb+=NegativeBinomialRandom(k,nmean);
1341 }
bb442091 1342 return (1-beta)*mulnp/2+beta*mulnb;
1343}
bb442091 1344
14244e14 1345//______________________________________________________________________________
7288f660 1346Double_t AliGlauberMC::GetdNdEtaTwoNBD ( const Double_t* p ) const
bb442091 1347{
1348 //Get particle density per unit of rapidity
1349 //using random numbers from two negative binomial distributions
7288f660 1350 //Parameters: k1 = related to distribition width of distribution 1=3
1351 // nmean1 = mean of distribution 1=4
1352 // k2 = related to distribition width of distribution 2=2
1353 // nmean2 = mean of distribution 2=11
1354 // alpha = set contributions of distrubitin 1 / distribution 2=0.4
1355 // beta = set contribution of participants / binary collisions to multiplicity =0.13
1356 Int_t k1 = TMath::Nint(p[0]);
1357 Double_t nmean1 = p[1];
1358 Int_t k2 = TMath::Nint(p[2]);
1359 Double_t nmean2 = p[3];
1360 Double_t alpha = p[4];
1361 Double_t beta = p[6];
1362
bb442091 1363 Double_t mulnp=0.0;
1364 for(Int_t i = 0; i<fNpart; i++)
1365 {
1366 mulnp+=DoubleNegativeBinomialRandom(k1,nmean1,k2,nmean2,alpha);
1367 }
1368 Double_t mulnb=0.0;
1369 for(Int_t i = 0; i<fNcoll; i++)
1370 {
1371 mulnb+=DoubleNegativeBinomialRandom(k1,nmean1,k2,nmean2,alpha);
1372 }
1373 Double_t mul=(1-beta)*mulnp/2+beta*mulnb;
1374 return mul;
bb442091 1375}
14244e14 1376
1b4934d5 1377//______________________________________________________________________________
1378Double_t AliGlauberMC::GetEccentricityPart() const
1379{
bb442091 1380 //get participant eccentricity of participants
a0278ee7 1381 if (fNpart<2) return 0.0;
2c0ff1e6 1382 return (TMath::Sqrt((fSy2Parts-fSx2Parts)*(fSy2Parts-fSx2Parts)+4*fSxyParts*fSxyParts)/(fSy2Parts+fSx2Parts));
1b4934d5 1383}
1384
106d1ca1 1385//_____________________________________________________________________________
bb442091 1386Double_t AliGlauberMC::GetEccentricityPartColl() const
1387{
1388 //get participant eccentricity of binary collisions
998de9b9 1389 if (fNcoll<0) return 0.0;
bb442091 1390 return (TMath::Sqrt((fSy2Coll-fSx2Coll)*(fSy2Coll-fSx2Coll)+4*fSxyColl*fSxyColl)/(fSy2Coll+fSx2Coll));
1391}
1392
4d264dfa 1393//_____________________________________________________________________________
1394Double_t AliGlauberMC::GetEccentricityPartCom() const
1395{
1396 //get participant eccentricity of binary collisions
998de9b9 1397 if (fNcom<0) return 0.0;
4d264dfa 1398 return (TMath::Sqrt((fSy2Com-fSx2Com)*(fSy2Com-fSx2Com)+4*fSxyCom*fSxyCom)/(fSy2Com+fSx2Com));
1399}
bb442091 1400//______________________________________________________________________________
1b4934d5 1401Double_t AliGlauberMC::GetEccentricity() const
1402{
bb442091 1403 //get standard eccentricity of participants
74d0740e 1404 if (fNpart<2) return 0.0;
2c0ff1e6 1405 return ((fSy2Parts-fSx2Parts)/(fSy2Parts+fSx2Parts));
1406}
1407
1408//______________________________________________________________________________
4d264dfa 1409Double_t AliGlauberMC::GetEccentricityColl() const
2c0ff1e6 1410{
4d264dfa 1411 //get standard eccentricity of binary collisions
998de9b9 1412 if (fNcoll<0) return 0.0;
4d264dfa 1413 return ((fSy2Coll-fSx2Coll)/(fSy2Coll+fSx2Coll));
1b4934d5 1414}
bb442091 1415
1416//______________________________________________________________________________
4d264dfa 1417Double_t AliGlauberMC::GetEccentricityCom() const
bb442091 1418{
ddf04a47 1419 //get standard eccentricity of combined weight
998de9b9 1420 if (fNcom<0) return 0.0;
1421 return ((fSy2Com-fSx2Com)/(fSy2Com+fSx2Com));
1422 //return (((fMeanY2Com-(fMeanYCom*fMeanYCom)) - (fMeanX2Com-(fMeanXCom*fMeanXCom))) / ((fMeanY2Com-(fMeanYCom*fMeanYCom)) + (fMeanX2Com-(fMeanXCom*fMeanXCom))) );
4d264dfa 1423}
1424
1425//______________________________________________________________________________
1426Double_t AliGlauberMC::GetStoa() const
1427{
1428 //get standard Transverse Overlap Area
1429 if (fNpart<2) return 0.0;
1430 return ( TMath::Pi()*(TMath::Sqrt(fSx2Parts))*(TMath::Sqrt(fSy2Parts)));
bb442091 1431}
14244e14 1432
1b4934d5 1433//______________________________________________________________________________
ec852657 1434Bool_t AliGlauberMC::NextEvent(Double_t bgen)
1435{
1b4934d5 1436 //Make a new event
bb442091 1437 Int_t nAttempts = 10; // set indices, max attempts and max comparisons
1438 Bool_t succes = kFALSE;
1439 for(Int_t j=0; j<nAttempts; j++)
1440 {
1441 if(bgen<0||!succes) //get impactparameter
1442 {
1443 bgen = TMath::Sqrt((fBMax*fBMax-fBMin*fBMin)*gRandom->Rndm()+fBMin*fBMin);
1444 }
97d98e48 1445 if ( (succes=CalcEvent(bgen)) ) break; //ends if we have particparts
bb442091 1446 }
1447 return succes;
ec852657 1448}
ec852657 1449//______________________________________________________________________________
106d1ca1 1450Double_t AliGlauberMC::GetEpsilon2Part() const
1451{
1452 //get participant eccentricity of participants
1453 if (fNpart<2) return 0.0;
1454 return (TMath::Sqrt(fMeanr2Cos2Phi*fMeanr2Cos2Phi+fMeanr2Sin2Phi*fMeanr2Sin2Phi)/fMeanr2);
1455}
106d1ca1 1456//______________________________________________________________________________
106d1ca1 1457Double_t AliGlauberMC::GetEpsilon3Part() const
1458{
1459 //get participant eccentricity of participants
1460 if (fNpart<2) return 0.0;
998de9b9 1461 return (TMath::Sqrt(fMeanr2Cos3Phi*fMeanr2Cos3Phi+fMeanr2Sin3Phi*fMeanr2Sin3Phi)/fMeanr2);
1462 //return (TMath::Sqrt(fMeanr3Cos3Phi*fMeanr3Cos3Phi+fMeanr3Sin3Phi*fMeanr3Sin3Phi)/fMeanr3);
106d1ca1 1463}
106d1ca1 1464//______________________________________________________________________________
106d1ca1 1465Double_t AliGlauberMC::GetEpsilon4Part() const
1466{
1467 //get participant eccentricity of participants
1468 if (fNpart<2) return 0.0;
998de9b9 1469 return (TMath::Sqrt(fMeanr2Cos4Phi*fMeanr2Cos4Phi+fMeanr2Sin4Phi*fMeanr2Sin4Phi)/fMeanr2);
1470 //return (TMath::Sqrt(fMeanr4Cos4Phi*fMeanr4Cos4Phi+fMeanr4Sin4Phi*fMeanr4Sin4Phi)/fMeanr4);
106d1ca1 1471}
1472
1473//______________________________________________________________________________
106d1ca1 1474Double_t AliGlauberMC::GetEpsilon5Part() const
1475{
1476 //get participant eccentricity of participants
1477 if (fNpart<2) return 0.0;
998de9b9 1478 return (TMath::Sqrt(fMeanr2Cos5Phi*fMeanr2Cos5Phi+fMeanr2Sin5Phi*fMeanr2Sin5Phi)/fMeanr2);
1479 //return (TMath::Sqrt(fMeanr5Cos5Phi*fMeanr5Cos5Phi+fMeanr5Sin5Phi*fMeanr5Sin5Phi)/fMeanr5);
106d1ca1 1480}
9a5780fa 1481
106d1ca1 1482//______________________________________________________________________________
9a5780fa 1483Double_t AliGlauberMC::GetEpsilon2Coll() const
1484{
1485 //get epsilon2 of binary collisions
998de9b9 1486 if (fNcoll<0) return 0.0;
9a5780fa 1487 return (TMath::Sqrt(fMeanr2Cos2PhiColl*fMeanr2Cos2PhiColl+fMeanr2Sin2PhiColl*fMeanr2Sin2PhiColl)/fMeanr2Coll);
1488}
1489
1490//______________________________________________________________________________
1491Double_t AliGlauberMC::GetEpsilon3Coll() const
1492{
7288f660 1493 //get epsilon3 of binary collisions
998de9b9 1494 if (fNcoll<0) return 0.0;
1495 return (TMath::Sqrt(fMeanr2Cos3PhiColl*fMeanr2Cos3PhiColl+fMeanr2Sin3PhiColl*fMeanr2Sin3PhiColl)/fMeanr2Coll);
1496 //return (TMath::Sqrt(fMeanr3Cos3PhiColl*fMeanr3Cos3PhiColl+fMeanr3Sin3PhiColl*fMeanr3Sin3PhiColl)/fMeanr3Coll);
9a5780fa 1497}
1498
1499//______________________________________________________________________________
1500Double_t AliGlauberMC::GetEpsilon4Coll() const
1501{
7288f660 1502 //get epsilon4 of binary collisions
998de9b9 1503 if (fNcoll<0) return 0.0;
1504 return (TMath::Sqrt(fMeanr2Cos4PhiColl*fMeanr2Cos4PhiColl+fMeanr2Sin4PhiColl*fMeanr2Sin4PhiColl)/fMeanr2Coll);
1505 //return (TMath::Sqrt(fMeanr4Cos4PhiColl*fMeanr4Cos4PhiColl+fMeanr4Sin4PhiColl*fMeanr4Sin4PhiColl)/fMeanr4Coll);
9a5780fa 1506}
106d1ca1 1507
9a5780fa 1508//______________________________________________________________________________
1509Double_t AliGlauberMC::GetEpsilon5Coll() const
1510{
1511 //get epsilon5 of binary collisions
998de9b9 1512 if (fNcoll<0) return 0.0;
1513 return (TMath::Sqrt(fMeanr2Cos5PhiColl*fMeanr2Cos5PhiColl+fMeanr2Sin5PhiColl*fMeanr2Sin5PhiColl)/fMeanr2Coll);
1514 //return (TMath::Sqrt(fMeanr5Cos5PhiColl*fMeanr5Cos5PhiColl+fMeanr5Sin5PhiColl*fMeanr5Sin5PhiColl)/fMeanr5Coll);
4d264dfa 1515}
1516
1517//______________________________________________________________________________
1518Double_t AliGlauberMC::GetEpsilon2Com() const
1519{
1520 //get epsilon2 of binary collisions
998de9b9 1521 if (fNcom<0) return 0.0;
4d264dfa 1522 return (TMath::Sqrt(fMeanr2Cos2PhiCom*fMeanr2Cos2PhiCom+fMeanr2Sin2PhiCom*fMeanr2Sin2PhiCom)/fMeanr2Com);
1523}
1524
1525//______________________________________________________________________________
1526Double_t AliGlauberMC::GetEpsilon3Com() const
1527{
1528 //get epsilon3 of binary collisions
998de9b9 1529 if (fNcom<0) return 0.0;
1530 return (TMath::Sqrt(fMeanr2Cos3PhiCom*fMeanr2Cos3PhiCom+fMeanr2Sin3PhiCom*fMeanr2Sin3PhiCom)/fMeanr2Com);
1531 //return (TMath::Sqrt(fMeanr3Cos3PhiCom*fMeanr3Cos3PhiCom+fMeanr3Sin3PhiCom*fMeanr3Sin3PhiCom)/fMeanr3Com);
4d264dfa 1532}
1533
1534//______________________________________________________________________________
1535Double_t AliGlauberMC::GetEpsilon4Com() const
1536{
1537 //get epsilon4 of binary collisions
998de9b9 1538 if (fNcom<0) return 0.0;
1539 return (TMath::Sqrt(fMeanr2Cos4PhiCom*fMeanr2Cos4PhiCom+fMeanr2Sin4PhiCom*fMeanr2Sin4PhiCom)/fMeanr2Com);
1540 //return (TMath::Sqrt(fMeanr4Cos4PhiCom*fMeanr4Cos4PhiCom+fMeanr4Sin4PhiCom*fMeanr4Sin4PhiCom)/fMeanr4Com);
4d264dfa 1541}
1542
1543//______________________________________________________________________________
1544Double_t AliGlauberMC::GetEpsilon5Com() const
1545{
1546 //get epsilon5 of binary collisions
998de9b9 1547 if (fNcom<0) return 0.0;
1548 return (TMath::Sqrt(fMeanr2Cos5PhiCom*fMeanr2Cos5PhiCom+fMeanr2Sin5PhiCom*fMeanr2Sin5PhiCom)/fMeanr2Com);
1549 //return (TMath::Sqrt(fMeanr5Cos5PhiCom*fMeanr5Cos5PhiCom+fMeanr5Sin5PhiCom*fMeanr5Sin5PhiCom)/fMeanr5Com);
9a5780fa 1550}
1551
1552//______________________________________________________________________________
8f9ee6e7 1553Double_t AliGlauberMC::GetPsi2() const
1554{
d10259b9 1555 return ((TMath::ATan2(fMeanr2Sin2Phi,fMeanr2Cos2Phi)+TMath::Pi())/2);
8f9ee6e7 1556}
1557
1558//______________________________________________________________________________
1559Double_t AliGlauberMC::GetPsi3() const
1560{
998de9b9 1561 return ((TMath::ATan2(fMeanr2Sin3Phi,fMeanr2Cos3Phi)+TMath::Pi())/3);
1562 //return ((TMath::ATan2(fMeanr3Sin3Phi,fMeanr3Cos3Phi)+TMath::Pi())/3);
8f9ee6e7 1563}
1564
1565//______________________________________________________________________________
1566Double_t AliGlauberMC::GetPsi4() const
1567{
998de9b9 1568 return ((TMath::ATan2(fMeanr2Sin4Phi,fMeanr2Cos4Phi)+TMath::Pi())/4);
1569 //return ((TMath::ATan2(fMeanr4Sin4Phi,fMeanr4Cos4Phi)+TMath::Pi())/4);
8f9ee6e7 1570}
1571
1572//______________________________________________________________________________
1573Double_t AliGlauberMC::GetPsi5() const
1574{
998de9b9 1575 return ((TMath::ATan2(fMeanr2Sin5Phi,fMeanr2Cos5Phi)+TMath::Pi())/5);
1576 //return ((TMath::ATan2(fMeanr5Sin5Phi,fMeanr5Cos5Phi)+TMath::Pi())/5);
8f9ee6e7 1577}
1578
4d264dfa 1579/*_____________________________________________________________________________
1580Double_t AliGlauberMC::GetE43Part() const
1581{
1582 //get participant eccentricity of participants
1583 if (fNpart<2) return 0.0;
1584 //return ((TMath::Sqrt(fMeanr2Cos2Phi*fMeanr2Cos2Phi+fMeanr2Sin2Phi*fMeanr2Sin2Phi)/fMeanr2) * (TMath::Sqrt(fMeanr2Cos2Phi*fMeanr2Cos2Phi+fMeanr2Sin2Phi*fMeanr2Sin2Phi)/fMeanr2) * (TMath::Sqrt(fMeanr2Cos4Phi*fMeanr2Cos4Phi+fMeanr2Sin4Phi*fMeanr2Sin4Phi)/fMeanr2) * TMath::Cos (((TMath::ATan2(fMeanr2Sin2Phi,fMeanr2Cos2Phi)+TMath::Pi())/2) - ((TMath::ATan2(fMeanr2Sin4Phi,fMeanr2Cos4Phi)+TMath::Pi())/4)));
1585 return ((TMath::Sqrt(fMeanr2Cos2Phi*fMeanr2Cos2Phi+fMeanr2Sin2Phi*fMeanr2Sin2Phi)/fMeanr2) * (TMath::Sqrt(fMeanr2Cos2Phi*fMeanr2Cos2Phi+fMeanr2Sin2Phi*fMeanr2Sin2Phi)/fMeanr2) * (TMath::Sqrt(fMeanr4Cos4Phi*fMeanr4Cos4Phi+fMeanr4Sin4Phi*fMeanr4Sin4Phi)/fMeanr4) * TMath::Cos(((TMath::ATan2(fMeanr2Sin2Phi,fMeanr2Cos2Phi)+TMath::Pi())/2) - ((TMath::ATan2(fMeanr4Sin4Phi,fMeanr4Cos4Phi)+TMath::Pi())/4)));
1586}
1587
1588
1589//______________________________________________________________________________
1590Double_t AliGlauberMC::GetE43Coll() const
1591{
1592 //get epsilon5 of binary collisions
1593 if (fNcoll<2) return 0.0;
1594 //return ((TMath::Sqrt(fMeanr2Cos2PhiColl*fMeanr2Cos2PhiColl+fMeanr2Sin2PhiColl*fMeanr2Sin2PhiColl)/fMeanr2Coll) * (TMath::Sqrt(fMeanr2Cos2PhiColl*fMeanr2Cos2PhiColl+fMeanr2Sin2PhiColl*fMeanr2Sin2PhiColl)/fMeanr2Coll) * (TMath::Sqrt(fMeanr2Cos4PhiColl*fMeanr2Cos4PhiColl+fMeanr2Sin4PhiColl*fMeanr2Sin4PhiColl)/fMeanr2Coll) * TMath::Cos(((TMath::ATan2(fMeanr2Sin2Phi,fMeanr2Cos2Phi)+TMath::Pi())/2) - ((TMath::ATan2(fMeanr2Sin4Phi,fMeanr2Cos4Phi)+TMath::Pi())/4)));
1595 return ((TMath::Sqrt(fMeanr2Cos2PhiColl*fMeanr2Cos2PhiColl+fMeanr2Sin2PhiColl*fMeanr2Sin2PhiColl)/fMeanr2Coll) * (TMath::Sqrt(fMeanr2Cos2PhiColl*fMeanr2Cos2PhiColl+fMeanr2Sin2PhiColl*fMeanr2Sin2PhiColl)/fMeanr2Coll) * (TMath::Sqrt(fMeanr4Cos4PhiColl*fMeanr4Cos4PhiColl+fMeanr4Sin4PhiColl*fMeanr4Sin4PhiColl)/fMeanr4Coll) * TMath::Cos(((TMath::ATan2(fMeanr2Sin2Phi,fMeanr2Cos2Phi)+TMath::Pi())/2) - ((TMath::ATan2(fMeanr4Sin4Phi,fMeanr4Cos4Phi)+TMath::Pi())/4)));
1596}
1597
1598//______________________________________________________________________________
1599Double_t AliGlauberMC::GetE43Com() const
1600{
1601 //get epsilon5 of binary collisions
1602 if (fNcom<2) return 0.0;
1603 //return ((TMath::Sqrt(fMeanr2Cos2PhiCom*fMeanr2Cos2PhiCom+fMeanr2Sin2PhiCom*fMeanr2Sin2PhiCom)/fMeanr2Com) * (TMath::Sqrt(fMeanr2Cos2PhiCom*fMeanr2Cos2PhiCom+fMeanr2Sin2PhiCom*fMeanr2Sin2PhiCom)/fMeanr2Com) * (TMath::Sqrt(fMeanr2Cos4PhiCom*fMeanr2Cos4PhiCom+fMeanr2Sin4PhiCom*fMeanr2Sin4PhiCom)/fMeanr2Com) * TMath::Cos(((TMath::ATan2(fMeanr2Sin2Phi,fMeanr2Cos2Phi)+TMath::Pi())/2) - ((TMath::ATan2(fMeanr2Sin4Phi,fMeanr2Cos4Phi)+TMath::Pi())/4)));
1604 return ((TMath::Sqrt(fMeanr2Cos2PhiCom*fMeanr2Cos2PhiCom+fMeanr2Sin2PhiCom*fMeanr2Sin2PhiCom)/fMeanr2Com) * (TMath::Sqrt(fMeanr2Cos2PhiCom*fMeanr2Cos2PhiCom+fMeanr2Sin2PhiCom*fMeanr2Sin2PhiCom)/fMeanr2Com) * (TMath::Sqrt(fMeanr4Cos4PhiCom*fMeanr4Cos4PhiCom+fMeanr4Sin4PhiCom*fMeanr4Sin4PhiCom)/fMeanr4Com) * TMath::Cos(((TMath::ATan2(fMeanr2Sin2Phi,fMeanr2Cos2Phi)+TMath::Pi())/2) - ((TMath::ATan2(fMeanr4Sin4Phi,fMeanr4Cos4Phi)+TMath::Pi())/4)));
1605}
1606//___________________________________________________________________________
1607
1608Double_t AliGlauberMC::GetPsi4m2() const
1609{
1610 //return (TMath::Cos(4*((TMath::ATan2(fMeanr2Sin4Phi,fMeanr2Cos4Phi)+TMath::Pi())/4) - ((TMath::ATan2(fMeanr2Sin2Phi,fMeanr2Cos2Phi)+TMath::Pi())/2)));
1611 return (TMath::Cos(4*(((TMath::ATan2(fMeanr4Sin4Phi,fMeanr4Cos4Phi)+TMath::Pi())/4)-((TMath::ATan2(fMeanr2Sin2Phi,fMeanr2Cos2Phi)+TMath::Pi())/2))));
1612}
1613*/
8f9ee6e7 1614//______________________________________________________________________________
ec852657 1615void AliGlauberMC::Run(Int_t nevents)
1616{
43215add 1617 //example run
bb442091 1618 cout << "Generating " << nevents << " events..." << endl;
1619 TString name(Form("nt_%s_%s",fANucleus.GetName(),fBNucleus.GetName()));
1620 TString title(Form("%s + %s (x-sect = %d mb)",fANucleus.GetName(),fBNucleus.GetName(),(Int_t) fXSect));
1621 if (fnt == 0)
1622 {
1623 fnt = new TNtuple(name,title,
ddf04a47 1624 "Npart:Ncoll:B:MeanX:MeanY:MeanX2:MeanY2:MeanXY:VarX:VarY:VarXY:MeanXSystem:MeanYSystem:MeanXA:MeanYA:MeanXB:MeanYB:VarE:Stoa:VarEColl:VarECom:VarEPart:VarEPartColl:VarEPartCom:dNdEta:dNdEtaGBW:dNdEtaTwoNBD:xsect:tAA:Epsl2:Epsl3:Epsl4:Epsl5:E2Coll:E3Coll:E4Coll:E5Coll:E2Com:E3Com:E4Com:E5Com:Psi2:Psi3:Psi4:Psi5");
bb442091 1625 fnt->SetDirectory(0);
1626 }
1627 Int_t q = 0;
1628 Int_t u = 0;
1629 for (Int_t i = 0; i<nevents; i++)
7288f660 1630 {
1631
1632 if(!NextEvent())
1633 {
1634 u++;
1635 continue;
1636 }
1637
1638 q++;
1639 //Float_t v[27];
4d264dfa 1640 Float_t v[45];
7288f660 1641 v[0] = GetNpart();
1642 v[1] = GetNcoll();
1643 v[2] = fBMC;
1644 v[3] = fMeanXParts;
1645 v[4] = fMeanYParts;
2c0ff1e6 1646 v[5] = fMeanX2Parts;
1647 v[6] = fMeanY2Parts;
1648 v[7] = fMeanXYParts;
1649 v[8] = fSx2Parts;
1650 v[9] = fSy2Parts;
1651 v[10] = fSxyParts;
7288f660 1652 v[11] = fMeanXSystem;
1653 v[12] = fMeanYSystem;
1654 v[13] = fMeanXA;
1655 v[14] = fMeanYA;
1656 v[15] = fMeanXB;
1657 v[16] = fMeanYB;
1658 v[17] = GetEccentricity();
2c0ff1e6 1659 v[18] = GetStoa();
1660 v[19] = GetEccentricityColl();
4d264dfa 1661 v[20] = GetEccentricityCom();
1662 v[21] = GetEccentricityPart();
1663 v[22] = GetEccentricityPartColl();
1664 v[23] = GetEccentricityPartCom();
7288f660 1665 if (fDoPartProd)
1666 {
4d264dfa 1667 v[24] = GetdNdEta();
1668 v[25] = GetdNdEta();
1669 v[26] = v[24]+v[25];
7288f660 1670 }
1671 else
bb442091 1672 {
2c0ff1e6 1673 v[24] = 0;
4d264dfa 1674 v[25] = 0;
1675 v[26] = 0;
106d1ca1 1676 }
4d264dfa 1677 v[27]=fXSect;
7288f660 1678
1679 Float_t mytAA=-999;
1680 if (GetNcoll()>0) mytAA=GetNcoll()/fXSect;
4d264dfa 1681 v[28]=mytAA;
7288f660 1682 //_____________epsilon2,3,4,4_______
4d264dfa 1683 v[29] = GetEpsilon2Part();
1684 v[30] = GetEpsilon3Part();
1685 v[31] = GetEpsilon4Part();
1686 v[32] = GetEpsilon5Part();
1687 v[33] = GetEpsilon2Coll();
1688 v[34] = GetEpsilon3Coll();
1689 v[35] = GetEpsilon4Coll();
1690 v[36] = GetEpsilon5Coll();
1691 v[37] = GetEpsilon2Com();
1692 v[38] = GetEpsilon3Com();
1693 v[39] = GetEpsilon4Com();
1694 v[40] = GetEpsilon5Com();
1695 v[41] = GetPsi2();
1696 v[42] = GetPsi3();
1697 v[43] = GetPsi4();
1698 v[44] = GetPsi5();
7288f660 1699 //always at the end
1700 fnt->Fill(v);
1701
1702 if ((i%100)==0) std::cout << "Generating Event # " << i << "... \r" << flush;
1703 }
bb442091 1704 std::cout << "Generating Event # " << nevents << "... \r" << endl << "Done! Succesfull events: " << q << " discarded events: " << u <<"."<< endl;
ec852657 1705}
1706
1707//---------------------------------------------------------------------------------
43215add 1708void AliGlauberMC::RunAndSaveNtuple( Int_t n,
1709 const Option_t *sysA,
1710 const Option_t *sysB,
ec852657 1711 Double_t signn,
1712 Double_t mind,
7288f660 1713 Double_t r,
1714 Double_t a,
ec852657 1715 const char *fname)
1716{
43215add 1717 //example run
566fc3d0 1718 AliGlauberMC mcg(sysA,sysB,signn);
1719 mcg.SetMinDistance(mind);
14244e14 1720 mcg.Setr(r);
1721 mcg.Seta(a);
566fc3d0 1722 mcg.Run(n);
1723 TNtuple *nt=mcg.GetNtuple();
bb442091 1724 TFile out(fname,"recreate",fname,9);
1725 if(nt) nt->Write();
14244e14 1726 printf("total cross section with a nucleon-nucleon cross section \t%f is \t%f",signn,mcg.GetTotXSect());
bb442091 1727 out.Close();
ec852657 1728}
1729
1730//---------------------------------------------------------------------------------
43215add 1731void AliGlauberMC::RunAndSaveNucleons( Int_t n,
1732 const Option_t *sysA,
bb442091 1733 Option_t *sysB,
43215add 1734 const Double_t signn,
ec852657 1735 Double_t mind,
1736 Bool_t verbose,
1737 const char *fname)
1738{
43215add 1739 //example run
566fc3d0 1740 AliGlauberMC mcg(sysA,sysB,signn);
1741 mcg.SetMinDistance(mind);
bb442091 1742 TFile *out=0;
1743 if(fname)
1744 out=new TFile(fname,"recreate",fname,9);
1745
1746 for(Int_t ievent=0; ievent<n; ievent++)
7288f660 1747 {
1748
1749 //get an event with at least one collision
1750 while(!mcg.NextEvent()) {}
1751
1752 //access, save and (if wanted) print out nucleons
1753 TObjArray* nucleons=mcg.GetNucleons();
1754 if(!nucleons) continue;
1755 if(out)
1756 nucleons->Write(Form("nucleonarray%d",ievent),TObject::kSingleKey);
1757
1758 if(verbose)
bb442091 1759 {
7288f660 1760 cout<<endl<<endl<<"EVENT NO: "<<ievent<<endl;
1761 cout<<"B = "<<mcg.GetB()<<" Npart = "<<mcg.GetNpart()<<endl<<endl;
1762 printf("Nucleus\t X\t Y\t Z\tNcoll\n");
1763 Int_t nNucls=nucleons->GetEntries();
1764 for(Int_t iNucl=0; iNucl<nNucls; iNucl++)
1765 {
1766 AliGlauberNucleon *nucl=(AliGlauberNucleon *)nucleons->UncheckedAt(iNucl);
1767 Char_t nucleus='A';
1768 if(nucl->IsInNucleusB()) nucleus='B';
1769 Double_t x=nucl->GetX();
1770 Double_t y=nucl->GetY();
1771 Double_t z=nucl->GetZ();
1772 Int_t ncoll=nucl->GetNColl();
1773 printf(" %c\t%2.2f\t%2.2f\t%2.2f\t%3d\n",nucleus,x,y,z,ncoll);
1774 }
bb442091 1775 }
7288f660 1776 }
bb442091 1777 if(out) delete out;
ec852657 1778}
1779
566fc3d0 1780//---------------------------------------------------------------------------------
1781void AliGlauberMC::Reset()
1782{
1783 //delete the ntuple
7288f660 1784 delete fnt;
1785 fnt=NULL;
566fc3d0 1786}