- Adapted comments for Doxygen
[u/mrichter/AliRoot.git] / THbtp / AliGenHBTprocessor.cxx
CommitLineData
36b81802 1/* Copyright(c) 1998-1999, ALICE Experiment at CERN, All rights reserved. *
2 * See cxx source for full Copyright notice */
3// Implementation of the interface for THBTprocessor
4// Author: Piotr Krzysztof Skowronski <Piotr.Skowronski@cern.ch>
5//////////////////////////////////////////////////////////////////////////////////
6//Wrapper class for "hbt processor" after burner
7//The origibal code is written in fortran by Lanny Ray
8//and is put in the directory $ALICE_ROOT/HBTprocessor
9//Detailed description is on the top of the file hbt_event_processor.f
10//
11//This class can be used ONLY with AliGenCocktailAfterBurner wrapper generator
12//i.e. (in Config.C)
13// ....
14// AliGenCocktailAfterBurner = gener = new AliGenCocktailAfterBurner();
15// gener->SetPhiRange(0, 360); //Set global parameters
16// gener->SetThetaRange(35., 145.);
17// AliGenHIJINGpara *hijing = new AliGenHIJINGpara(10000); //Add some generator
18// hijing->SetMomentumRange(0, 999);
19// gener->AddGenerator(hijing,"HIJING PARAMETRIZATION",1);
20//
21// AliGenHBTprocessor *hbtp = new AliGenHBTprocessor(); //create object
22// hbtp->SetRefControl(2); //set parameters
23// hbtp->SetSwitch1D(1);
24// hbtp->SetR0(6);//fm - radius
25// hbtp->SetLambda(0.7);//chaocity parameter
26// gener->AddAfterBurner(hbtp,"HBT PROCESSOR",1); //add to the main generator
27//
28// gener->Init();
29//
30//CAUTIONS:
31// A) HBT PROCESSOR NEEDS MORE THAN ONE EVENT TO WORK
32// AS MORE AS IT BETTER WORKS
33// B) IT IS ABLE TO "ADD" CORRELATIONS ONLY UP TO TWO PARTICLE TYPES AT ONES
88cb7938 34//
35//
36// Artificial particle dennsity enhancment feature
37// HBT Processor is unable to process correctly low multiplicity particles
38// even if the high statistics (more events) is supplied.
39// For that reason was implemented artificial multiplicity enhancement
40// - see also comments in HBT Processor source file (HBTP/hbt_event_processor.f)
41// or web page (http://alisoft.cern.ch/people/skowron/hbtprocessor/hbteventprocessor.html)
42// concerning trk_accep or SetTrackRejectionFactor method in this class
43// Fortran is cheated by masking several events a single one. Number is defined by fEventMerge
44// variable. In case it is equal to 1, there is no masking and the feature is off.
45// But, for example if it is 5, and we have 100 events, number of events passed to fortran is 100/5=20
46// When fortran asks about multiplcity of event, let sey, 1 - multiplicity of events 5 to 9 is summed up
47// and returned.
48//
49//
36b81802 50//////////////////////////////////////////////////////////////////////////////////
51
52// 11.11.2001 Piotr Skowronski
53// Setting seed (date) in RNG in the constructor
54
55// 09.10.2001 Piotr Skowronski
56//
57// Theta - Eta cohernecy facilities added:
58// AliGenerator::SetThetaRange method overriden
59// Static methods added
60// EtaToTheta
61// ThetaToEta
62// DegreesToRadians
63// RadiansToDegrees
64//
65// Class description comments put on proper place
66
67// 27.09.2001 Piotr Skowronski
68// removing of redefinition of defaults velues
69// in method's implementation.
70//
71//
72
73#include "AliGenHBTprocessor.h"
88cb7938 74
88cb7938 75#include <TParticle.h>
76#include "THBTprocessor.h"
77
36b81802 78#include "AliStack.h"
36b81802 79#include "AliGenCocktailAfterBurner.h"
371660fe 80#include "AliLog.h"
36b81802 81
82
83
84ClassImp(AliGenHBTprocessor)
85
88cb7938 86Int_t AliGenHBTprocessor::fgDebug = 0;
b456eb2e 87static TRandom* gAliRandom;//RNG to be used by the fortran code
36b81802 88
89AliGenCocktailAfterBurner* GetGenerator();
90/*******************************************************************/
88cb7938 91AliGenHBTprocessor::AliGenHBTprocessor(const AliGenHBTprocessor& in):
b456eb2e 92 AliGenerator(in),
93 fHBTprocessor(0x0),
94 fHbtPStatCodes(0x0),
95 fEventMerge(1)
36b81802 96{
97//copy contructor
98 // AliGenHBTprocessor::AliGenHBTprocessor();
36b81802 99
b456eb2e 100}
101AliGenHBTprocessor::AliGenHBTprocessor():
102 AliGenerator(),
103 fHBTprocessor(0x0),
104 fHbtPStatCodes(0x0),
105 fEventMerge(1)
36b81802 106{
107 //
108 // Standard constructor
109 // Sets default veues of all parameters
36b81802 110
111 SetName("AliGenHBTprocessor");
112 SetTitle("AliGenHBTprocessor");
113
b456eb2e 114 gAliRandom = fRandom;
36b81802 115 fHBTprocessor = new THBTprocessor();
116
117 fNPDGCodes = 0;
118 DefineParticles();
119
120 SetTrackRejectionFactor();
121 SetRefControl();
122 SetPIDs();
123 SetNPIDtypes();
124 SetDeltap();
125 SetMaxIterations();
126 SetDelChi();
127 SetIRand();
128 SetLambda();
129 SetR1d() ;
130 SetRSide();
131 SetROut() ;
132 SetRLong() ;
133 SetRPerp();
134 SetRParallel();
135 SetR0();
136 SetQ0();
137 SetSwitch1D();
138 SetSwitch3D();
139 SetSwitchType();
140 SetSwitchCoherence();
141 SetSwitchCoulomb();
142 SetSwitchFermiBose();
143 //SetMomentumRange();
144 SetPtRange();
145 SetPxRange();
146 SetPyRange();
147 SetPzRange();
148 SetPhiRange();
149 SetEtaRange();
150 SetNPtBins();
151 SetNPhiBins();
152 SetNEtaBins();
153 SetNPxBins();
154 SetNPyBins();
155 SetNPzBins();
156 SetNBins1DFineMesh();
157 SetBinSize1DFineMesh();
158 SetNBins1DCoarseMesh();
159 SetBinSize1DCoarseMesh();
160 SetNBins3DFineMesh();
161 SetBinSize3DFineMesh();
162 SetNBins3DCoarseMesh();
163 SetBinSize3DCoarseMesh();
164 SetNBins3DFineProjectMesh();
165}
166
167/*******************************************************************/
168
169
170/*******************************************************************/
171
172AliGenHBTprocessor::~AliGenHBTprocessor()
173{
174//destructor
175 CleanStatusCodes();
176 if (fHBTprocessor) delete fHBTprocessor; //delete generator
177
178}
179
180/*******************************************************************/
181
182void AliGenHBTprocessor::InitStatusCodes()
183{
184 //creates and inits status codes array to zero
185 AliGenCocktailAfterBurner *cab = GetGenerator();
186
187 if(!cab) Fatal("InitStatusCodes()","Can not find AliGenCocktailAfterBurner generator");
188
189 Int_t nev = cab->GetNumberOfEvents();
190
191 fHbtPStatCodes = new Int_t* [nev];
192 for( Int_t i =0; i<nev;i++)
193 {
194 Int_t nprim = cab->GetStack(i)->GetNprimary();
195 fHbtPStatCodes[i] = new Int_t[nprim];
196 for (Int_t k =0 ;k<nprim;k++)
197 fHbtPStatCodes[i][k] =0;
198
199 }
200
201}
202/*******************************************************************/
203
204void AliGenHBTprocessor::CleanStatusCodes()
205{
206 //Cleans up status codes
207 if (fHbtPStatCodes)
208 {
209 for (Int_t i =0; i<GetGenerator()->GetNumberOfEvents(); i++)
210 delete [] fHbtPStatCodes[i];
211 delete fHbtPStatCodes;
212 fHbtPStatCodes = 0;
213 }
214
215}
216/*******************************************************************/
217
218void AliGenHBTprocessor::Init()
219 {
220 //sets up parameters in generator
221
222 THBTprocessor *thbtp = fHBTprocessor;
223
224
225 thbtp->SetTrackRejectionFactor(fTrackRejectionFactor);
226 thbtp->SetRefControl(fReferenceControl);
227
228 if ((fPid[0] == fPid[1]) || (fPid[0] == 0) || (fPid[1] == 0))
229 {
230 if (fPid[0] == 0)
231 thbtp->SetPIDs(IdFromPDG(fPid[1]) ,0);
232 else
233 thbtp->SetPIDs(IdFromPDG(fPid[0]) ,0);
234 thbtp->SetNPIDtypes(1);
235
236 if (fSwitchType !=1)
237 Warning("AliGenHBTprocessor::Init","\nThere is only one particle type set,\n\
238 and Switch_Type differnt then 1,\n which does not make sense.\n\
239 Setting it to 1.\n");
240
241 SetSwitchType(1);
242 }
243 else
244 {
245 thbtp->SetPIDs(IdFromPDG(fPid[0]) ,IdFromPDG(fPid[1]));
246 SetNPIDtypes(2);
247 thbtp->SetSwitchType(fSwitchType);
248 }
249
250 thbtp->SetMaxIterations(fMaxit);
251 thbtp->SetDelChi(fDelchi);
252 thbtp->SetIRand(fIrand);
253 thbtp->SetLambda(fLambda);
254 thbtp->SetR1d(fR1d);
255 thbtp->SetRSide(fRside);
256 thbtp->SetROut(fRout);
257 thbtp->SetRLong(fRlong);
258 thbtp->SetRPerp(fRperp);
259 thbtp->SetRParallel(fRparallel);
260 thbtp->SetR0(fR0);
261 thbtp->SetQ0(fQ0);
262 thbtp->SetSwitch1D(fSwitch1d);
263 thbtp->SetSwitch3D(fSwitch3d);
264 thbtp->SetSwitchType(fSwitchType);
265 thbtp->SetSwitchCoherence(fSwitchCoherence);
266 thbtp->SetSwitchCoulomb(fSwitchCoulomb);
267 thbtp->SetSwitchFermiBose(fSwitchFermiBose);
268 thbtp->SetPtRange(fPtMin,fPtMax);
269 thbtp->SetPxRange(fPxMin,fPxMax);
270 thbtp->SetPyRange(fPyMin,fPyMax);
271 thbtp->SetPzRange(fPzMin,fPzMax);
88cb7938 272 thbtp->SetPhiRange(fPhiMin*180./((Float_t)TMath::Pi())+180.0, //casting is because if fPhiMin = 180.0 then
273 fPhiMax*180./((Float_t)TMath::Pi())+180.0);//TMath::Pi() != TMath::Pi()*fPhiMin/180.0,
36b81802 274 thbtp->SetEtaRange(fEtaMin,fEtaMax);
275 thbtp->SetNPtBins(fNPtBins);
276 thbtp->SetNPhiBins(fNPhiBins);
277 thbtp->SetNEtaBins(fNEtaBins);
278 thbtp->SetNPxBins(fNPxBins);
279 thbtp->SetNPyBins(fNPyBins);
280 thbtp->SetNPzBins(fNPzBins);
281 thbtp->SetNBins1DFineMesh(fN1dFine);
282 thbtp->SetBinSize1DFineMesh(fBinsize1dFine);
283 thbtp->SetNBins1DCoarseMesh(fN1dCoarse);
284 thbtp->SetBinSize1DCoarseMesh(fBinsize1dCoarse);
285 thbtp->SetNBins3DFineMesh(fN3dFine);
286 thbtp->SetBinSize3DFineMesh(fBinsize3dFine);
287 thbtp->SetNBins3DCoarseMesh(fN3dCoarse);
288 thbtp->SetBinSize3DCoarseMesh(fBinsize3dCoarse);
289 thbtp->SetNBins3DFineProjectMesh(fN3dFineProject);
88cb7938 290
291 thbtp->SetPrintFull(fPrintFull);
36b81802 292
293 }
294/*******************************************************************/
295
296void AliGenHBTprocessor::Generate()
297 {
298 //starts processig
299 AliGenCocktailAfterBurner* cab = GetGenerator();
300 if (cab == 0x0)
301 {
302 Fatal("Generate()","AliGenHBTprocessor needs AliGenCocktailAfterBurner to be main generator");
303 }
304 if (cab->GetNumberOfEvents() <2)
305 {
306 Fatal("Generate()",
307 "HBT Processor needs more than 1 event to work properly,\
308 but there is only %d set", cab->GetNumberOfEvents());
309 }
310
311
312 fHBTprocessor->Initialize(); //reset all fields of common blocks
313 //in case there are many HBT processors
314 //run one after one (i.e. in AliCocktailAfterBurner)
315 //between Init() called and Generate there might
316 Init(); //be different instance running - be sure that we have our settings
317
318 InitStatusCodes(); //Init status codes
319
320 fHBTprocessor->GenerateEvent(); //Generates event
321
322 CleanStatusCodes(); //Clean Status codes - thet are not needed anymore
323 }
324
325/*******************************************************************/
326
327
328/*******************************************************************/
329void AliGenHBTprocessor::GetParticles(TClonesArray * particles) const
330 {
331 //practically dumm
332 if (!particles)
333 {
88cb7938 334 Error("GetParticles","Particles has to be initialized");
36b81802 335 return;
336 }
337 fHBTprocessor->ImportParticles(particles);
338 }
339
340/*******************************************************************/
341
342Int_t AliGenHBTprocessor::GetHbtPStatusCode(Int_t part) const
343{
344//returns the status code of the given particle in the active event
345//see SetActiveEvent in the bottom of AliGenHBTprocessor.cxx
346//and in AliCocktailAfterBurner
88cb7938 347 Int_t ev, idx;
348 GetTrackEventIndex(part,ev,idx);
349 if ( (ev<0) || (idx<0) )
350 {
351 Error("GetHbtPStatusCode","GetTrackEventIndex returned error");
352 return 0;
353 }
0798b21e 354 return fHbtPStatCodes[ev][idx];
88cb7938 355
36b81802 356}
357
358/*******************************************************************/
359void AliGenHBTprocessor::SetHbtPStatusCode(Int_t hbtstatcode, Int_t part)
360{
361 //Sets the given status code to given particle number (part) in the active event
88cb7938 362 Int_t ev, idx;
363 GetTrackEventIndex(part,ev,idx);
364 if ( (ev<0) || (idx<0) )
365 {
366 Error("GetHbtPStatusCode","GetTrackEventIndex returned error");
367 return;
368 }
369 else fHbtPStatCodes[ev][idx] = hbtstatcode;
36b81802 370}
371
372/*******************************************************************/
373
374void AliGenHBTprocessor::SetTrackRejectionFactor(Float_t trf) //def 1.0
375 {
376 //Sets the Track Rejection Factor
377 //variates in range 0.0 <-> 1.0
378 //Describes the factor of particles rejected from the output.
379 //Used only in case of low muliplicity particles e.g. lambdas.
380 //Processor generates addisional particles and builds the
381 //correletions on such a statistics.
382 //At the end these particels are left in the event according
383 //to this factor: 1==all particles are left
384 // 0==all are removed
385
386 fTrackRejectionFactor=trf;
387 fHBTprocessor->SetTrackRejectionFactor(trf);
388 }
389/*******************************************************************/
390
391void AliGenHBTprocessor::SetRefControl(Int_t rc) //default 2
392 {
393 //Sets the Refernce Control Switch
394 //switch wether read reference histograms from file =1
395 // compute from input events =2 - default
396 fReferenceControl=rc;
397 fHBTprocessor->SetRefControl(rc);
398 }
399/*******************************************************************/
400
401void AliGenHBTprocessor::SetPIDs(Int_t pid1,Int_t pid2)
402 {
403 //default pi+ pi-
404 //Sets PDG Codes of particles to be processed, default \\Pi^{+} and \\Pi^{-}
405 //This method accepts PDG codes which is
406 //in opposite to THBBProcessor which accepts GEANT PID
407 if ( (pid1 == 0) && (pid2 == 0) )
408 {
409 Error("AliGenHBTprocessor::SetPIDs","Sensless Particle Types setting: 0 0, Ignoring\n");
410 }
411
412 fPid[0]=pid1;
413 fPid[1]=pid2;
414
415 if(pid1 == pid2)
416 {
417 fHBTprocessor->SetPIDs(IdFromPDG(pid1) ,0);
418 SetNPIDtypes(1);
419 SetSwitchType(1);
420 }
421 else
422 {
423 fHBTprocessor->SetPIDs(IdFromPDG(pid1) ,IdFromPDG(pid2));
424 SetNPIDtypes(2);
425 }
426 }
427/*******************************************************************/
428
429void AliGenHBTprocessor::SetNPIDtypes(Int_t npidt)
430 {
431 //Number ofparticle types to be processed - default 2
432 //see AliGenHBTprocessor::SetPIDs
433
434 fNPidTypes = npidt;
435 fHBTprocessor->SetNPIDtypes(npidt);
436 }
437/*******************************************************************/
438
439void AliGenHBTprocessor::SetDeltap(Float_t deltp)
440 {
441 //default = 0.1 GeV
442 //maximum range for random momentum shifts in GeV/c;
443 //px,py,pz independent; Default = 0.1 GeV/c.
444 fDeltap=deltp;
445 fHBTprocessor->SetDeltap(deltp);
446 }
447/*******************************************************************/
448
449void AliGenHBTprocessor::SetMaxIterations(Int_t maxiter)
450 {
451 //maximum # allowed iterations thru all the
452 //tracks for each event. Default = 50.
453 //If maxit=0, then calculate the correlations
454 //for the input set of events without doing the
455 //track adjustment procedure.
456
457 fMaxit=maxiter;
458 fHBTprocessor->SetMaxIterations(maxiter);
459 }
460
461/*******************************************************************/
462void AliGenHBTprocessor::SetDelChi(Float_t dc)
463 {
464 //min % change in total chi-square for which
465 //the track adjustment iterations may stop,
466 //Default = 0.1%.
467
468 fDelchi=dc;
469 fHBTprocessor->SetDelChi(dc);
470 }
471/*******************************************************************/
472
473void AliGenHBTprocessor::SetIRand(Int_t irnd)
474 {
475 //if fact dummy - used only for compatibility
476 //we are using AliRoot built in RNG
477 //dummy in fact since we are using aliroot build-in RNG
478 //Sets renaodom number generator
479 fIrand=irnd;
480 fHBTprocessor->SetIRand(irnd);
481 }
482/*******************************************************************/
483
484void AliGenHBTprocessor::SetLambda(Float_t lam)
485 {
486 //default = 0.6
487 // Sets Chaoticity Parameter
488 fLambda = lam;
489 fHBTprocessor->SetLambda(lam);
490 }
491/*******************************************************************/
492void AliGenHBTprocessor::SetR1d(Float_t r)
493 {
494 //default 7.0
495 //Sets Spherical source model radius (fm)
496 fR1d = r;
497 fHBTprocessor->SetR1d(r);
498 }
499/*******************************************************************/
500void AliGenHBTprocessor::SetRSide(Float_t rs)
501 {
502 //default rs = 6.0
503 //Rside,Rout,Rlong = Non-spherical Bertsch-Pratt source model (fm)
504
505 fRside = rs;
506 fHBTprocessor->SetRSide(rs);
507 }
508/*******************************************************************/
509void AliGenHBTprocessor::SetROut(Float_t ro)
510 {
511 //default ro = 7.0
512 //Rside,Rout,Rlong = Non-spherical Bertsch-Pratt source model (fm)
513 fRout = ro;
514 fHBTprocessor->SetROut(ro);
515 }
516/*******************************************************************/
517void AliGenHBTprocessor::SetRLong(Float_t rl)
518 {
519 //default rl = 4.0
520 //Rside,Rout,Rlong = Non-spherical Bertsch-Pratt source model (fm)
521 fRlong = rl;
522 fHBTprocessor->SetRLong(rl);
523 }
524/*******************************************************************/
525void AliGenHBTprocessor::SetRPerp(Float_t rp)
526 {
527 //default rp = 6.0
528 //Rperp,Rparallel,R0= Non-spherical Yano-Koonin-Podgoretski source model (fm).
529 fRperp = rp;
530 fHBTprocessor->SetRPerp(rp);
531 }
532/*******************************************************************/
533void AliGenHBTprocessor::SetRParallel(Float_t rprl)
534 {
535 //default rprl = 4.0
536 //Rperp,Rparallel,R0= Non-spherical Yano-Koonin-Podgoretski source model (fm).
537 fRparallel = rprl;
538 fHBTprocessor->SetRParallel(rprl);
539 }
540/*******************************************************************/
541void AliGenHBTprocessor::SetR0(Float_t r0)
542 {
543 //default r0 = 4.0
544 //Rperp,Rparallel,R0= Non-spherical Yano-Koonin-Podgoretski source model (fm).
545 fR0 = r0;
546 fHBTprocessor->SetR0(r0);
547 }
548/*******************************************************************/
549void AliGenHBTprocessor::SetQ0(Float_t q0)
550 {
551 //default q0 = 9.0
552 //Sets Q0 = NA35 Coulomb parameter for finite source size in (GeV/c)
553 // if fSwitchCoulomb = 2
554 // = Spherical Coulomb source radius in (fm)
555 // if switchCoulomb = 3, used to interpolate the
556 // input Pratt/Cramer discrete Coulomb source
557 // radii tables.
558 fQ0 = q0;
559 fHBTprocessor->SetQ0(q0);
560 }
561
562/*******************************************************************/
563void AliGenHBTprocessor::SetSwitch1D(Int_t s1d)
564 {
565//default s1d = 3
566// Sets fSwitch1d
567// = 0 to not compute the 1D two-body //orrelations.
568// = 1 to compute this using Q-invariant
569// = 2 to compute this using Q-total
570// = 3 to compute this using Q-3-ve//tor
571
572 fSwitch1d = s1d;
573 fHBTprocessor->SetSwitch1D(s1d);
574 }
575/*******************************************************************/
576void AliGenHBTprocessor::SetSwitch3D(Int_t s3d)
577 {
578//default s3d = 0
579// Sets fSwitch3d
580// = 0 to not compute the 3D two-body correlations.
581// = 1 to compute this using the side-out-long form
582// = 2 to compute this using the Yanno-Koonin-Pogoredskij form.
583
584 fSwitch3d = s3d;
585 fHBTprocessor->SetSwitch3D(s3d);
586 }
587/*******************************************************************/
588void AliGenHBTprocessor::SetSwitchType(Int_t st)
589 {
590//default st = 3
591// Sets switch_type = 1 to fit only the like pair correlations
592// = 2 to fit only the unlike pair correlations
593// = 3 to fit both the like and unlike pair correl.
594//See SetPIDs and Init
595//If only one particle type is set, unly==1 makes sens
596
597 fSwitchType = st;
598 fHBTprocessor->SetSwitchType(st);
599 }
600/*******************************************************************/
601void AliGenHBTprocessor::SetSwitchCoherence(Int_t sc)
602 {
603// default sc = 0
604// switchCoherence = 0 for purely incoherent source (but can have
605// lambda < 1.0)
606// = 1 for mixed incoherent and coherent source
607
608 fSwitchCoherence = sc;
609 fHBTprocessor->SetSwitchCoherence(sc);
610 }
611/*******************************************************************/
612void AliGenHBTprocessor::SetSwitchCoulomb(Int_t scol)
613 {
614//default scol = 2
615// switchCoulomb = 0 no Coulomb correction
616// = 1 Point source Gamow correction only
617// = 2 NA35 finite source size correction
618// = 3 Pratt/Cramer finite source size correction;
619// interpolated from input tables.
620 fSwitchCoulomb =scol;
621 fHBTprocessor->SetSwitchCoulomb(scol);
622 }
623/*******************************************************************/
624void AliGenHBTprocessor::SetSwitchFermiBose(Int_t sfb)
625 {
626//default sfb = 1
627// switchFermiBose = 1 Boson pairs
628// = -1 Fermion pairs
629
630 fSwitchFermiBose = sfb;
631 fHBTprocessor->SetSwitchFermiBose(sfb);
632 }
633/*******************************************************************/
634void AliGenHBTprocessor::SetPtRange(Float_t ptmin, Float_t ptmax)
635 {
636// default ptmin = 0.1, ptmax = 0.98
637//Sets Pt range (GeV)
638 AliGenerator::SetPtRange(ptmin,ptmax);
639 fHBTprocessor->SetPtRange(ptmin,ptmax);
640 }
641
642/*******************************************************************/
643void AliGenHBTprocessor::SetPxRange(Float_t pxmin, Float_t pxmax)
644 {
645//default pxmin = -1.0, pxmax = 1.0
646//Sets Px range
647 fPxMin =pxmin;
648 fPxMax =pxmax;
649 fHBTprocessor->SetPxRange(pxmin,pxmax);
650 }
651/*******************************************************************/
652void AliGenHBTprocessor::SetPyRange(Float_t pymin, Float_t pymax)
653 {
654//default pymin = -1.0, pymax = 1.0
655//Sets Py range
656 fPyMin =pymin;
657 fPyMax =pymax;
658 fHBTprocessor->SetPyRange(pymin,pymax);
659 }
660/*******************************************************************/
661void AliGenHBTprocessor::SetPzRange(Float_t pzmin, Float_t pzmax)
662 {
663//default pzmin = -3.6, pzmax = 3.6
664//Sets Py range
665 fPzMin =pzmin;
666 fPzMax =pzmax;
667 fHBTprocessor->SetPzRange(pzmin,pzmax);
668 }
371660fe 669void AliGenHBTprocessor::SetMomentumRange(Float_t /*pmin*/, Float_t /*pmax*/)
36b81802 670 {
671 //default pmin=0, pmax=0
672 //Do not use this method!
673 MayNotUse("AliGenHBTprocessor::SetMomentumRange Method is Dummy");
674 }
675
676 /*******************************************************************/
677void AliGenHBTprocessor::SetPhiRange(Float_t phimin, Float_t phimax)
678 {
679//
680//Sets \\Phi range
681 AliGenerator::SetPhiRange(phimin,phimax);
682
88cb7938 683 fHBTprocessor->SetPhiRange(phimin+180.0,phimax+180.0);
36b81802 684 }
685/*******************************************************************/
686void AliGenHBTprocessor::SetEtaRange(Float_t etamin, Float_t etamax)
687 {
688//default etamin = -1.5, etamax = 1.5
689//Sets \\Eta range
690 fEtaMin= etamin;
691 fEtaMax =etamax;
692 fHBTprocessor->SetEtaRange(etamin,etamax);
693
694 //set the azimothal angle range in the AliGeneraor -
695 //to keep coherency between azimuthal angle and pseudorapidity
696 //DO NOT CALL this->SetThetaRange, because it calls this method (where we are)
697 //which must cause INFINITE LOOP
698 AliGenerator::SetThetaRange(RadiansToDegrees(EtaToTheta(fEtaMin)),
699 RadiansToDegrees(EtaToTheta(fEtaMax)));
700
701 }
702/*******************************************************************/
703void AliGenHBTprocessor::SetThetaRange(Float_t thetamin, Float_t thetamax)
704{
705 //default thetamin = 0, thetamax = 180
706 //Azimuthal angle, override AliGenerator method which uses widely (i.e. wrapper generators)
707 //core fortran HBTProcessor uses Eta (pseudorapidity)
708 //so these methods has to be synchronized
709
710 AliGenerator::SetThetaRange(thetamin,thetamax);
711
712 SetEtaRange( ThetaToEta(fThetaMin) , ThetaToEta(fThetaMax) );
713
714}
715
716/*******************************************************************/
717void AliGenHBTprocessor::SetNPtBins(Int_t nptbin)
718 {
719 //default nptbin = 50
720 //set number of Pt bins
721 fNPtBins= nptbin;
722 fHBTprocessor->SetNPtBins(nptbin);
723 }
724/*******************************************************************/
725void AliGenHBTprocessor::SetNPhiBins(Int_t nphibin)
726{
727 //default nphibin = 50
728 //set number of Phi bins
729 fNPhiBins=nphibin;
730 fHBTprocessor->SetNPhiBins(nphibin);
731}
732/*******************************************************************/
733void AliGenHBTprocessor::SetNEtaBins(Int_t netabin)
734{
735 //default netabin = 50
736 //set number of Eta bins
737 fNEtaBins = netabin;
738 fHBTprocessor->SetNEtaBins(netabin);
739}
740/*******************************************************************/
741void AliGenHBTprocessor::SetNPxBins(Int_t npxbin)
742{
743 //default npxbin = 20
744 //set number of Px bins
745 fNPxBins = npxbin;
746 fHBTprocessor->SetNPxBins(npxbin);
747}
748/*******************************************************************/
749void AliGenHBTprocessor::SetNPyBins(Int_t npybin)
750{
751 //default npybin = 20
752 //set number of Py bins
753 fNPyBins = npybin;
754 fHBTprocessor->SetNPyBins(npybin);
755}
756/*******************************************************************/
757void AliGenHBTprocessor::SetNPzBins(Int_t npzbin)
758{
759 //default npzbin = 70
760 //set number of Pz bins
761 fNPzBins = npzbin;
762 fHBTprocessor->SetNPzBins(npzbin);
763}
764/*******************************************************************/
765void AliGenHBTprocessor::SetNBins1DFineMesh(Int_t n)
766{
767//default n = 10
768//Sets the number of bins in the 1D mesh
769 fN1dFine =n;
770 fHBTprocessor->SetNBins1DFineMesh(n);
771
772}
773/*******************************************************************/
774void AliGenHBTprocessor::SetBinSize1DFineMesh(Float_t x)
775{
776//default x=0.01
777//Sets the bin size in the 1D mesh
778 fBinsize1dFine = x;
779 fHBTprocessor->SetBinSize1DFineMesh(x);
780}
781/*******************************************************************/
782
783void AliGenHBTprocessor::SetNBins1DCoarseMesh(Int_t n)
784{
785//default n =2
786//Sets the number of bins in the coarse 1D mesh
787 fN1dCoarse =n;
788 fHBTprocessor->SetNBins1DCoarseMesh(n);
789}
790/*******************************************************************/
791void AliGenHBTprocessor::SetBinSize1DCoarseMesh(Float_t x)
792{
793//default x=0.05
794//Sets the bin size in the coarse 1D mesh
795 fBinsize1dCoarse =x;
796 fHBTprocessor->SetBinSize1DCoarseMesh(x);
797}
798/*******************************************************************/
799
800void AliGenHBTprocessor::SetNBins3DFineMesh(Int_t n)
801{
802//default n = 8
803//Sets the number of bins in the 3D mesh
804 fN3dFine =n;
805 fHBTprocessor->SetNBins3DFineMesh(n);
806}
807/*******************************************************************/
808void AliGenHBTprocessor::SetBinSize3DFineMesh(Float_t x)
809{
810//default x=0.01
811//Sets the bin size in the 3D mesh
812 fBinsize3dFine =x;
813 fHBTprocessor->SetBinSize3DFineMesh(x);
814}
815/*******************************************************************/
816
817void AliGenHBTprocessor::SetNBins3DCoarseMesh(Int_t n)
818{
819//default n = 2
820//Sets the number of bins in the coarse 3D mesh
821
822 fN3dCoarse = n;
823 fHBTprocessor->SetNBins3DCoarseMesh(n);
824}
825/*******************************************************************/
826void AliGenHBTprocessor::SetBinSize3DCoarseMesh(Float_t x)
827{
828//default x=0.08
829//Sets the bin size in the coarse 3D mesh
830 fBinsize3dCoarse = x;
831 fHBTprocessor->SetBinSize3DCoarseMesh(x);
832}
833/*******************************************************************/
834
835void AliGenHBTprocessor::SetNBins3DFineProjectMesh(Int_t n )
836{
837//default n =3
838//Sets the number of bins in the fine project mesh
839 fN3dFineProject = n;
840 fHBTprocessor->SetNBins3DFineProjectMesh(n);
841}
842/*******************************************************************/
88cb7938 843void AliGenHBTprocessor::SetPrintFull(Int_t flag)
844{
b456eb2e 845//sets the print full flag
88cb7938 846 fPrintFull = flag;
847 fHBTprocessor->SetPrintFull(flag);
848}
849
850
851/*******************************************************************/
852
853Int_t AliGenHBTprocessor::GetNumberOfEvents()
854{
b456eb2e 855//returns number of available events
88cb7938 856 AliGenerator* g = gAlice->Generator();
857 AliGenCocktailAfterBurner* cab = (g)?dynamic_cast<AliGenCocktailAfterBurner*>(g):0x0;
858 if (cab == 0x0)
859 {
860 Fatal("GetNumberOfEvents","Master Generator is not an AliGenCocktailAfterBurner");
861 return 0;
862 }
863 return (Int_t)TMath::Ceil(cab->GetNumberOfEvents()/((Float_t)fEventMerge));
864}
865
866/*******************************************************************/
36b81802 867
88cb7938 868void AliGenHBTprocessor::SetActiveEventNumber(Int_t n)
869{
b456eb2e 870//sets the active event
88cb7938 871 fActiveStack = n*fEventMerge;
371660fe 872 AliDebug(1,Form("Settimg active event %d passed %d",fActiveStack,n));
88cb7938 873}
874/*******************************************************************/
36b81802 875
88cb7938 876Int_t AliGenHBTprocessor::GetNumberOfTracks()
877{
b456eb2e 878//returns number of tracks in active event
88cb7938 879 AliGenerator* g = gAlice->Generator();
880 AliGenCocktailAfterBurner* cab = (g)?dynamic_cast<AliGenCocktailAfterBurner*>(g):0x0;
881 if (cab == 0x0)
882 {
883 Fatal("GetNumberOfEvents","Master Generator is not an AliGenCocktailAfterBurner");
884 return 0;
885 }
886 Int_t n = 0;
887 for (Int_t i = fActiveStack;i < fActiveStack+fEventMerge; i++)
888 {
889 if (i >= GetNumberOfEvents()) break; //protection not to overshoot nb of events
890 AliStack* stack = cab->GetStack(i);
891 if (stack == 0x0)
892 Error("GetNumberOfTracks","There is no stack %d",i);
893
894 n+=stack->GetNprimary();
895 }
896 return n;
897}
898/*******************************************************************/
899
900void AliGenHBTprocessor::SetNEventsToMerge(Int_t nev)
901{
b456eb2e 902 //sets number of events to merge
88cb7938 903 if (nev > 0 ) fEventMerge = nev;
904}
905/*******************************************************************/
906
907TParticle* AliGenHBTprocessor::GetTrack(Int_t n)
908{
909//returns track that hbtp thinks is n in active event
371660fe 910 AliDebug(5,Form("n = %d",n));
88cb7938 911 AliGenerator* g = gAlice->Generator();
912 AliGenCocktailAfterBurner* cab = (g)?dynamic_cast<AliGenCocktailAfterBurner*>(g):0x0;
913 if (cab == 0x0)
914 {
915 Fatal("GetTrackEventIndex","Master Generator is not an AliGenCocktailAfterBurner");
916 return 0;
917 }
918
919 Int_t ev, idx;
920 GetTrackEventIndex(n,ev,idx);
371660fe 921 AliDebug(5,Form("Event = %d Particle = %d",ev,idx));
88cb7938 922 if ( (ev<0) || (idx<0) )
923 {
924 Error("GetTrack","GetTrackEventIndex returned error");
925 return 0x0;
926 }
371660fe 927 AliDebug(5,Form("Number of Tracks in Event(%d) = %d",ev,cab->GetStack(ev)->GetNprimary()));
0798b21e 928 return cab->GetStack(ev)->Particle(idx);//safe - in case stack does not exist
88cb7938 929}
36b81802 930/*******************************************************************/
931
88cb7938 932void AliGenHBTprocessor::GetTrackEventIndex(Int_t n, Int_t &evno, Int_t &index) const
933{
934 //returns event(stack) number and particle index
935 AliGenerator* g = gAlice->Generator();
936 AliGenCocktailAfterBurner* cab = (g)?dynamic_cast<AliGenCocktailAfterBurner*>(g):0x0;
937 if (cab == 0x0)
938 {
939 Fatal("GetTrackEventIndex","Master Generator is not an AliGenCocktailAfterBurner");
940 return;
941 }
942
943 evno = -1;
944 index = -1;
945 for (Int_t i = fActiveStack;i < fActiveStack+fEventMerge; i++)
946 {
947 AliStack* stack = cab->GetStack(i);
948 if (stack == 0x0)
949 {
950 Error("GetTrackEventIndex","There is no stack %d",i);
951 return;
952 }
953
954 Int_t ntracks = stack->GetNprimary();
371660fe 955 AliDebug(10,Form("Event %d has %d tracks. n = %d",i,ntracks,n));
88cb7938 956
957 if ( ntracks > n)
958 {
959 evno = i;
960 index = n;
961 return ;
962 }
963 else
964 {
965 n-=ntracks;
966 continue;
967 }
968 }
969 Error("GetTrackEventIndex","Could not find given track");
970}
971
972/*******************************************************************/
973/*******************************************************************/
974/*******************************************************************/
36b81802 975
976
977
978
979
980void AliGenHBTprocessor::DefineParticles()
981{
982 //
983 // Load standard numbers for GEANT particles and PDG conversion
984 fNPDGCodes = 0; //this is done in the constructor - but in any case ...
985
986 fPDGCode[fNPDGCodes++]=-99; // 0 = unused location
987 fPDGCode[fNPDGCodes++]=22; // 1 = photon
988 fPDGCode[fNPDGCodes++]=-11; // 2 = positron
989 fPDGCode[fNPDGCodes++]=11; // 3 = electron
990 fPDGCode[fNPDGCodes++]=12; // 4 = neutrino e
991 fPDGCode[fNPDGCodes++]=-13; // 5 = muon +
992 fPDGCode[fNPDGCodes++]=13; // 6 = muon -
993 fPDGCode[fNPDGCodes++]=111; // 7 = pi0
994 fPDGCode[fNPDGCodes++]=211; // 8 = pi+
995 fPDGCode[fNPDGCodes++]=-211; // 9 = pi-
996 fPDGCode[fNPDGCodes++]=130; // 10 = Kaon Long
997 fPDGCode[fNPDGCodes++]=321; // 11 = Kaon +
998 fPDGCode[fNPDGCodes++]=-321; // 12 = Kaon -
999 fPDGCode[fNPDGCodes++]=2112; // 13 = Neutron
1000 fPDGCode[fNPDGCodes++]=2212; // 14 = Proton
1001 fPDGCode[fNPDGCodes++]=-2212; // 15 = Anti Proton
1002 fPDGCode[fNPDGCodes++]=310; // 16 = Kaon Short
1003 fPDGCode[fNPDGCodes++]=221; // 17 = Eta
1004 fPDGCode[fNPDGCodes++]=3122; // 18 = Lambda
1005 fPDGCode[fNPDGCodes++]=3222; // 19 = Sigma +
1006 fPDGCode[fNPDGCodes++]=3212; // 20 = Sigma 0
1007 fPDGCode[fNPDGCodes++]=3112; // 21 = Sigma -
1008 fPDGCode[fNPDGCodes++]=3322; // 22 = Xi0
1009 fPDGCode[fNPDGCodes++]=3312; // 23 = Xi-
1010 fPDGCode[fNPDGCodes++]=3334; // 24 = Omega-
1011 fPDGCode[fNPDGCodes++]=-2112; // 25 = Anti Proton
1012 fPDGCode[fNPDGCodes++]=-3122; // 26 = Anti Proton
1013 fPDGCode[fNPDGCodes++]=-3222; // 27 = Anti Sigma -
1014 fPDGCode[fNPDGCodes++]=-3212; // 28 = Anti Sigma 0
1015 fPDGCode[fNPDGCodes++]=-3112; // 29 = Anti Sigma 0
1016 fPDGCode[fNPDGCodes++]=-3322; // 30 = Anti Xi 0
1017 fPDGCode[fNPDGCodes++]=-3312; // 31 = Anti Xi +
1018 fPDGCode[fNPDGCodes++]=-3334; // 32 = Anti Omega +
1019}
1020
1021/*******************************************************************/
1022Int_t AliGenHBTprocessor::IdFromPDG(Int_t pdg) const
1023{
1024 //
1025 // Return Geant3 code from PDG and pseudo ENDF code
1026 //
1027 for(Int_t i=0;i<fNPDGCodes;++i)
1028 if(pdg==fPDGCode[i]) return i;
1029 return -1;
1030}
1031Int_t AliGenHBTprocessor::PDGFromId(Int_t id) const
1032{
1033 //
1034 // Return PDG code and pseudo ENDF code from Geant3 code
1035 //
1036 if(id>0 && id<fNPDGCodes) return fPDGCode[id];
1037 else return -1;
1038}
1039Double_t AliGenHBTprocessor::ThetaToEta(Double_t arg)
1040 {
1041 //converts etha(azimuthal angle) to Eta (pseudorapidity). Argument in radians
88cb7938 1042
1043 if(arg>= TMath::Pi()) return 708.0;//This number is the biggest wich not crashes exp(x)
1044 if(arg<= 0.0) return -708.0;//
1045 if(arg == TMath::Pi()/2.) return 0.0;//
36b81802 1046
36b81802 1047 if (arg > 0.0)
1048 {
88cb7938 1049 return TMath::Log( TMath::Tan(arg/2.)) ;
36b81802 1050 }
1051 else
88cb7938 1052 {
1053 return -TMath::Log( TMath::Tan(-arg/2.)) ;
36b81802 1054 }
1055 }
1056
1057/*******************************************************************/
1058/****** ROUTINES USED FOR COMMUNUCATION ********/
1059/******************** WITH FORTRAN ********************/
1060/*******************************************************************/
1061
1062#ifndef WIN32
1063 # define hbtpran hbtpran_
1064 # define alihbtp_puttrack alihbtp_puttrack_
1065 # define alihbtp_gettrack alihbtp_gettrack_
1066 # define alihbtp_getnumberevents alihbtp_getnumberevents_
1067 # define alihbtp_getnumbertracks alihbtp_getnumbertracks_
1068 # define alihbtp_initialize alihbtp_initialize_
1069 # define alihbtp_setactiveeventnumber alihbtp_setactiveeventnumber_
1070 # define alihbtp_setparameters alihbtp_setparameters_
1071 # define type_ofCall
1072
1073#else
1074 # define hbtpran HBTPRAN
1075 # define alihbtp_puttrack ALIHBTP_PUTTRACK
1076 # define alihbtp_gettrack ALIHBTP_GETTRACK
1077 # define alihbtp_getnumberevents ALIHBTP_GETNUMBEREVENTS
1078 # define alihbtp_getnumbertracks ALIHBTP_GETNUMBERTRACKS
1079 # define alihbtp_initialize ALIHBTP_INITIALIZE
1080 # define alihbtp_setactiveeventnumber ALIHBTP_SETACTIVEEVENTNUMBER
1081 # define alihbtp_setparameters ALIHBTP_SETPARAMETERS
1082 # define type_ofCall _stdcall
1083#endif
1084
36b81802 1085/*******************************************************************/
1086
1087AliGenCocktailAfterBurner* GetGenerator()
1088 {
1089 // This function has two tasks:
1090 // Check if environment is OK (exist gAlice and generator)
1091 // Returns pointer to genrator
1092 //to be changed with TFolders
1093
1094 if(!gAlice)
1095 {
88cb7938 1096 ::Error("AliGenHBTprocessor.cxx: GetGenerator()",
1097 "There is NO gALICE! Check what you are doing!");
1098 ::Fatal("AliGenHBTprocessor.cxx: GetGenerator()",
1099 "Running HBT Processor without gAlice... Exiting \n");
1100 return 0x0;//pro forma
36b81802 1101 }
1102 AliGenerator * gen = gAlice->Generator();
1103
1104 if (!gen)
1105 {
88cb7938 1106 ::Fatal("AliGenHBTprocessor.cxx: GetGenerator()",
1107 "There is no generator in gAlice, exiting\n");
36b81802 1108 return 0x0;
1109 }
1110
1111 //we do not sure actual type of the genetator
1112 //and simple casting is risky - we use ROOT machinery to do safe cast
1113
1114 TClass* cabclass = AliGenCocktailAfterBurner::Class(); //get AliGenCocktailAfterBurner TClass
1115 TClass* genclass = gen->IsA();//get TClass of the generator we got from galice
1116 //use casting implemented in TClass
1117 //cast gen to cabclass
1118 AliGenCocktailAfterBurner* cab=(AliGenCocktailAfterBurner*)genclass->DynamicCast(cabclass,gen);
1119
1120 if (cab == 0x0)//if generator that we got is not AliGenCocktailAfterBurner or its descendant we get null
1121 { //then quit with error
88cb7938 1122 ::Fatal("AliGenHBTprocessor.cxx: GetGenerator()",
1123 "The main Generator is not a AliGenCocktailAfterBurner, exiting\n");
36b81802 1124 return 0x0;
1125 }
36b81802 1126 return cab;
36b81802 1127 }
1128/*******************************************************************/
1129
1130AliGenHBTprocessor* GetAliGenHBTprocessor()
1131{
1132//returns pointer to the current instance of AliGenHBTprocessor in
1133//AliGenCocktailAfterBurner (can be more than one)
1134//
1135 AliGenCocktailAfterBurner* gen = GetGenerator();
1136 AliGenerator* g = gen->GetCurrentGenerator();
1137 if(g == 0x0)
1138 {
88cb7938 1139 ::Fatal("AliGenHBTprocessor.cxx: GetAliGenHBTprocessor()",
36b81802 1140 "Can not get the current generator. Exiting");
88cb7938 1141 return 0x0;//pro forma
36b81802 1142 }
1143
1144 TClass* hbtpclass = AliGenHBTprocessor::Class(); //get AliGenCocktailAfterBurner TClass
1145 TClass* gclass = g->IsA();//get TClass of the current generator we got from CAB
1146 AliGenHBTprocessor* hbtp = (AliGenHBTprocessor*)gclass->DynamicCast(hbtpclass,g);//try to cast
1147 if (hbtp == 0x0)
1148 {
88cb7938 1149 ::Fatal("AliGenHBTprocessor.cxx: GetAliGenHBTprocessor()",
1150 "\nCurrernt generator in AliGenCocktailAfterBurner is not a AliGenHBTprocessor, exiting\n");
36b81802 1151 return 0x0;
1152 }
1153 return hbtp;
1154}
1155
1156/*******************************************************************/
1157extern "C" void type_ofCall alihbtp_setparameters()
1158 {
1159 //dummy
1160 }
1161
1162extern "C" void type_ofCall alihbtp_initialize()
1163 {
1164 //dummy
1165 }
1166
1167/*******************************************************************/
1168
1169extern "C" void type_ofCall alihbtp_getnumberevents(Int_t &nev)
1170 {
371660fe 1171 //passes number of events to the fortran
b456eb2e 1172 if(AliGenHBTprocessor::GetDebug())
371660fe 1173 AliInfoGeneral("alihbtp_getnumberevents",Form("(%d) ....",nev));
88cb7938 1174
1175 AliGenHBTprocessor* gen = GetAliGenHBTprocessor();//we dont check because it is done in function
1176 nev = gen->GetNumberOfEvents();
36b81802 1177
b456eb2e 1178 if(AliGenHBTprocessor::GetDebug()>5)
371660fe 1179 AliInfoGeneral("alihbtp_getnumberevents",Form("EXITED N Ev = %d",nev));
36b81802 1180 }
1181
1182/*******************************************************************/
1183
1184extern "C" void type_ofCall alihbtp_setactiveeventnumber(Int_t & nev)
1185 {
1186//sets active event in generator (AliGenCocktailAfterBurner)
1187
b456eb2e 1188 if(AliGenHBTprocessor::GetDebug()>5)
88cb7938 1189 ::Info("AliGenHBTpocessor.cxx: alihbtp_setactiveeventnumber","(%d)",nev);
b456eb2e 1190 if(AliGenHBTprocessor::GetDebug()>0)
88cb7938 1191 ::Info("AliGenHBTpocessor.cxx: alihbtp_setactiveeventnumber","Asked for event %d",nev-1);
1192
1193 AliGenHBTprocessor* gen = GetAliGenHBTprocessor();//we dont check because it is done in function
36b81802 1194
88cb7938 1195 gen->SetActiveEventNumber(nev - 1); //fortran numerates events from 1 to N
36b81802 1196
b456eb2e 1197 if(AliGenHBTprocessor::GetDebug()>5)
88cb7938 1198 ::Info("AliGenHBTpocessor.cxx: alihbtp_setactiveeventnumber","EXITED returned %d",nev);
36b81802 1199 }
1200/*******************************************************************/
1201
1202extern "C" void type_ofCall alihbtp_getnumbertracks(Int_t &ntracks)
1203 {
1204//passes number of particles in active event to the fortran
b456eb2e 1205 if(AliGenHBTprocessor::GetDebug()>5)
88cb7938 1206 ::Info("AliGenHBTpocessor.cxx: alihbtp_getnumbertracks","(%d)",ntracks);
36b81802 1207
88cb7938 1208 AliGenHBTprocessor* gen = GetAliGenHBTprocessor();//we dont check because it is done in function
36b81802 1209
88cb7938 1210 ntracks = gen->GetNumberOfTracks();
b456eb2e 1211 if(AliGenHBTprocessor::GetDebug()>5)
88cb7938 1212 ::Info("AliGenHBTpocessor.cxx: alihbtp_getnumbertracks","EXITED Ntracks = %d",ntracks);
36b81802 1213 }
1214
1215/*******************************************************************/
1216
1217extern "C" void type_ofCall
1218 alihbtp_puttrack(Int_t & n,Int_t& flag, Float_t& px,
1219 Float_t& py, Float_t& pz, Int_t& geantpid)
1220 {
1221//sets new parameters (momenta) in track number n
1222// in the active event
1223// n - number of the track in active event
1224// flag - flag of the track
1225// px,py,pz - momenta
1226// geantpid - type of the particle - Geant Particle ID
1227
b456eb2e 1228 if(AliGenHBTprocessor::GetDebug()>5)
88cb7938 1229 ::Info("AliGenHBTpocessor.cxx: alihbtp_puttrack","(%d)",n);
36b81802 1230
88cb7938 1231 AliGenHBTprocessor* gen = GetAliGenHBTprocessor();//we dont check because it is done in function
36b81802 1232
88cb7938 1233 TParticle * track = gen->GetTrack(n-1);
1234 if (track == 0x0)
1235 {
1236 ::Error("AliGenHBTprocessor.cxx","Can not get track from AliGenHBTprocessor");
1237 return;
1238 }
36b81802 1239
1240 //check to be deleted
88cb7938 1241 if (geantpid != (gen->IdFromPDG( track->GetPdgCode() )))
36b81802 1242 {
88cb7938 1243 ::Error("AliGenHBTprocessor.cxx: alihbtp_puttrack",
1244 "SOMETHING IS GOING BAD:\n GEANTPIDS ARE NOT THE SAME");
36b81802 1245 }
1246
b456eb2e 1247 if(AliGenHBTprocessor::GetDebug()>0)
36b81802 1248 if (px != track->Px())
88cb7938 1249 ::Info("AliGenHBTprocessor.cxx: alihbtp_puttrack",
1250 "Px diff. = %f", px - track->Px());
36b81802 1251
b456eb2e 1252 if(AliGenHBTprocessor::GetDebug()>3)
88cb7938 1253 ::Info("AliGenHBTprocessor.cxx: alihbtp_puttrack",
1254 "track->GetPdgCode() --> %d",track->GetPdgCode());
36b81802 1255
1256
1257
1258 Float_t m = track->GetMass();
1259 Float_t e = TMath::Sqrt(m*m+px*px+py*py+pz*pz);
1260 track->SetMomentum(px,py,pz,e);
1261
88cb7938 1262 gen->SetHbtPStatusCode(flag,n-1);
36b81802 1263
b456eb2e 1264 if(AliGenHBTprocessor::GetDebug()>5) ::Info("AliGenHBTprocessor.cxx: alihbtp_puttrack","EXITED ");
36b81802 1265 }
1266
1267/*******************************************************************/
1268
1269extern "C" void type_ofCall
1270 alihbtp_gettrack(Int_t & n,Int_t & flag, Float_t & px,
1271 Float_t & py, Float_t & pz, Int_t & geantpid)
1272
1273 {
1274//passes track parameters to the fortran
1275// n - number of the track in active event
1276// flag - flag of the track
1277// px,py,pz - momenta
1278// geantpid - type of the particle - Geant Particle ID
1279
b456eb2e 1280 if(AliGenHBTprocessor::GetDebug()>5) ::Info("AliGenHBTprocessor.cxx: alihbtp_gettrack","(%d)",n);
36b81802 1281 AliGenHBTprocessor* g = GetAliGenHBTprocessor();
88cb7938 1282
1283 TParticle * track = g->GetTrack(n-1);
36b81802 1284
1285 flag = g->GetHbtPStatusCode(n-1);
1286
1287 px = (Float_t)track->Px();
1288 py = (Float_t)track->Py();
1289 pz = (Float_t)track->Pz();
1290
1291 geantpid = g->IdFromPDG( track->GetPdgCode() );
1292
b456eb2e 1293 if(AliGenHBTprocessor::GetDebug()>5) ::Info("AliGenHBTprocessor.cxx: alihbtp_gettrack","EXITED");
36b81802 1294 }
1295
1296/*******************************************************************/
1297extern "C" Float_t type_ofCall hbtpran(Int_t &)
1298{
1299//interface to the random number generator
b456eb2e 1300 return gAliRandom->Rndm();
36b81802 1301}
1302
1303/*******************************************************************/
1304
1305
1306/*******************************************************************/