]> git.uio.no Git - u/mrichter/AliRoot.git/blame - EVGEN/AliGenParam.cxx
git://git.uio.no / u / mrichter / AliRoot.git / blame
? search:
summary | shortlog | log | commit | commitdiff | tree
blob | blame (incremental) | history | HEAD
update from Asako's tasks
[u/mrichter/AliRoot.git] / EVGEN / AliGenParam.cxx
CommitLineData
4c039060 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
88cb7938 16/* $Id$ */
6d4dd317 17
18// Class to generate particles from using paramtrized pT and y distributions.
19// Distributions are obtained from pointer to object of type AliGenLib.
20// (For example AliGenMUONlib)
21// Decays are performed using Pythia.
22// andreas.morsch@cern.ch
23
ac3faee4 24#include <TCanvas.h>
7ca4655f 25#include <TClonesArray.h>
fa480368 26#include <TDatabasePDG.h>
ac3faee4 27#include <TF1.h>
28#include <TH1F.h>
fa480368 29#include <TLorentzVector.h>
ac3faee4 30#include <TMath.h>
31#include <TParticle.h>
32#include <TParticlePDG.h>
828fff0d 33#include <TROOT.h>
ac3faee4 34#include <TVirtualMC.h>
fa480368 35
ac3faee4 36#include "AliDecayer.h"
37#include "AliGenMUONlib.h"
38#include "AliGenParam.h"
5d12ce38 39#include "AliMC.h"
ac3faee4 40#include "AliRun.h"
7cf36cae 41#include "AliGenEventHeader.h"
fe4da5cc 42
43ClassImp(AliGenParam)
44
45//------------------------------------------------------------
46
6078e216 47//Begin_Html
48/*
1439f98e 49 <img src="picts/AliGenParam.gif">
6078e216 50*/
51//End_Html
fe4da5cc 52
fe4da5cc 53//____________________________________________________________
6078e216 54AliGenParam::AliGenParam()
55: fPtParaFunc(0),
1c56e311 56 fYParaFunc(0),
57 fIpParaFunc(0),
6078e216 58 fV2ParaFunc(0),
1c56e311 59 fPtPara(0),
60 fYPara(0),
6078e216 61 fV2Para(0),
62 fdNdPhi(0),
1c56e311 63 fParam(0),
64 fdNdy0(0.),
65 fYWgt(0.),
66 fPtWgt(0.),
67 fBias(0.),
68 fTrials(0),
69 fDeltaPt(0.01),
18e09c20 70 fSelectAll(kFALSE),
1c56e311 71 fDecayer(0)
fe4da5cc 72{
6078e216 73 // Default constructor
b22ee262 74}
2ad38d56 75//____________________________________________________________
f2c13e03 76AliGenParam::AliGenParam(Int_t npart, const AliGenLib * Library, Int_t param, const char* tname)
1c56e311 77 :AliGenMC(npart),
78 fPtParaFunc(Library->GetPt(param, tname)),
79 fYParaFunc (Library->GetY (param, tname)),
80 fIpParaFunc(Library->GetIp(param, tname)),
6078e216 81 fV2ParaFunc(Library->GetV2(param, tname)),
1c56e311 82 fPtPara(0),
83 fYPara(0),
6078e216 84 fV2Para(0),
85 fdNdPhi(0),
1c56e311 86 fParam(param),
87 fdNdy0(0.),
88 fYWgt(0.),
89 fPtWgt(0.),
90 fBias(0.),
91 fTrials(0),
92 fDeltaPt(0.01),
18e09c20 93 fSelectAll(kFALSE),
1c56e311 94 fDecayer(0)
b22ee262 95{
6078e216 96 // Constructor using number of particles parameterisation id and library
8b31bfac 97 fName = "Param";
98 fTitle= "Particle Generator using pT and y parameterisation";
34f60c01 99 fAnalog = kAnalog;
b22ee262 100 SetForceDecay();
fe4da5cc 101}
fe4da5cc 102//____________________________________________________________
1c56e311 103AliGenParam::AliGenParam(Int_t npart, Int_t param, const char* tname, const char* name):
104 AliGenMC(npart),
105 fPtParaFunc(0),
106 fYParaFunc (0),
107 fIpParaFunc(0),
6078e216 108 fV2ParaFunc(0),
1c56e311 109 fPtPara(0),
110 fYPara(0),
6078e216 111 fV2Para(0),
112 fdNdPhi(0),
1c56e311 113 fParam(param),
114 fdNdy0(0.),
115 fYWgt(0.),
116 fPtWgt(0.),
117 fBias(0.),
118 fTrials(0),
119 fDeltaPt(0.01),
18e09c20 120 fSelectAll(kFALSE),
1c56e311 121 fDecayer(0)
fe4da5cc 122{
6078e216 123 // Constructor using parameterisation id and number of particles
124 //
1c56e311 125 fName = name;
126 fTitle= "Particle Generator using pT and y parameterisation";
8b31bfac 127
675e9664 128 AliGenLib* pLibrary = new AliGenMUONlib();
675e9664 129 fPtParaFunc = pLibrary->GetPt(param, tname);
130 fYParaFunc = pLibrary->GetY (param, tname);
131 fIpParaFunc = pLibrary->GetIp(param, tname);
6078e216 132 fV2ParaFunc = pLibrary->GetV2(param, tname);
b22ee262 133
34f60c01 134 fAnalog = kAnalog;
b22ee262 135 fChildSelect.Set(5);
136 for (Int_t i=0; i<5; i++) fChildSelect[i]=0;
137 SetForceDecay();
138 SetCutOnChild();
139 SetChildMomentumRange();
140 SetChildPtRange();
141 SetChildPhiRange();
142 SetChildThetaRange();
fe4da5cc 143}
2ad38d56 144//____________________________________________________________
fe4da5cc 145
34f60c01 146AliGenParam::AliGenParam(Int_t npart, Int_t param,
75e0cc59 147 Double_t (*PtPara) (const Double_t*, const Double_t*),
148 Double_t (*YPara ) (const Double_t* ,const Double_t*),
6078e216 149 Double_t (*V2Para) (const Double_t* ,const Double_t*),
65fb704d 150 Int_t (*IpPara) (TRandom *))
1c56e311 151 :AliGenMC(npart),
152
153 fPtParaFunc(PtPara),
154 fYParaFunc(YPara),
155 fIpParaFunc(IpPara),
6078e216 156 fV2ParaFunc(V2Para),
1c56e311 157 fPtPara(0),
158 fYPara(0),
6078e216 159 fV2Para(0),
160 fdNdPhi(0),
1c56e311 161 fParam(param),
162 fdNdy0(0.),
163 fYWgt(0.),
164 fPtWgt(0.),
165 fBias(0.),
166 fTrials(0),
167 fDeltaPt(0.01),
18e09c20 168 fSelectAll(kFALSE),
1c56e311 169 fDecayer(0)
886b6f73 170{
6078e216 171 // Constructor
172 // Gines Martinez 1/10/99
8b31bfac 173 fName = "Param";
174 fTitle= "Particle Generator using pT and y parameterisation";
175
34f60c01 176 fAnalog = kAnalog;
886b6f73 177 fChildSelect.Set(5);
178 for (Int_t i=0; i<5; i++) fChildSelect[i]=0;
179 SetForceDecay();
180 SetCutOnChild();
749070b6 181 SetChildMomentumRange();
182 SetChildPtRange();
183 SetChildPhiRange();
184 SetChildThetaRange();
886b6f73 185}
186
fe4da5cc 187//____________________________________________________________
188AliGenParam::~AliGenParam()
189{
6078e216 190 // Destructor
fe4da5cc 191 delete fPtPara;
192 delete fYPara;
6078e216 193 delete fV2Para;
194 delete fdNdPhi;
fe4da5cc 195}
196
197//____________________________________________________________
198void AliGenParam::Init()
199{
6078e216 200 // Initialisation
fa480368 201
2904363f 202 if (gMC) fDecayer = gMC->GetDecayer();
fe4da5cc 203 //Begin_Html
204 /*
1439f98e 205 <img src="picts/AliGenParam.gif">
fe4da5cc 206 */
207 //End_Html
c5c3e4b0 208 char name[256];
8f8399ab 209 snprintf(name, 256, "pt-parameterisation for %s", GetName());
c5c3e4b0 210
2ad38d56 211 if (fPtPara) fPtPara->Delete();
c5c3e4b0 212 fPtPara = new TF1(name, fPtParaFunc, fPtMin, fPtMax,0);
828fff0d 213 gROOT->GetListOfFunctions()->Remove(fPtPara);
6078e216 214 // Set representation precision to 10 MeV
c5c3e4b0 215 Int_t npx= Int_t((fPtMax - fPtMin) / fDeltaPt);
749070b6 216
217 fPtPara->SetNpx(npx);
c5c3e4b0 218
8f8399ab 219 snprintf(name, 256, "y-parameterisation for %s", GetName());
2ad38d56 220 if (fYPara) fYPara->Delete();
c5c3e4b0 221 fYPara = new TF1(name, fYParaFunc, fYMin, fYMax, 0);
828fff0d 222 gROOT->GetListOfFunctions()->Remove(fYPara);
223
6078e216 224 snprintf(name, 256, "v2-parameterisation for %s", GetName());
225 if (fV2Para) fV2Para->Delete();
226 fV2Para = new TF1(name, fV2ParaFunc, fPtMin, fPtMax, 0);
227 // fV2Para = new TF1(name, "2*[0]/(1+TMath::Exp([1]*([2]-x)))-[0]", fPtMin, fPtMax);
228 // fV2Para->SetParameter(0, 0.236910);
229 // fV2Para->SetParameter(1, 1.71122);
230 // fV2Para->SetParameter(2, 0.0827617);
231 //gROOT->GetListOfFunctions()->Remove(fV2Para); //TR: necessary?
232
233 snprintf(name, 256, "dNdPhi for %s", GetName());
234 if (fdNdPhi) fdNdPhi->Delete();
235 fdNdPhi = new TF1(name, "1+2*[0]*TMath::Cos(2*(x-[1]))", fPhiMin, fPhiMax);
236 //gROOT->GetListOfFunctions()->Remove(fdNdPhi); //TR: necessary?
749070b6 237
8f8399ab 238 snprintf(name, 256, "pt-for-%s", GetName());
c5c3e4b0 239 TF1 ptPara(name ,fPtParaFunc, 0, 15, 0);
8f8399ab 240 snprintf(name, 256, "y-for-%s", GetName());
c5c3e4b0 241 TF1 yPara(name, fYParaFunc, -6, 6, 0);
b7601ac4 242
6078e216 243 //
244 // dN/dy| y=0
749070b6 245 Double_t y1=0;
246 Double_t y2=0;
247
248 fdNdy0=fYParaFunc(&y1,&y2);
6078e216 249 //
250 // Integral over generation region
faf00b0d 251#if ROOT_VERSION_CODE < ROOT_VERSION(5,99,0)
389a9124 252 Float_t intYS = yPara.Integral(fYMin, fYMax,(Double_t*) 0x0,1.e-6);
253 Float_t intPt0 = ptPara.Integral(0,15,(Double_t *) 0x0,1.e-6);
254 Float_t intPtS = ptPara.Integral(fPtMin,fPtMax,(Double_t*) 0x0,1.e-6);
faf00b0d 255#else
256 Float_t intYS = yPara.Integral(fYMin, fYMax,1.e-6);
257 Float_t intPt0 = ptPara.Integral(0,15,1.e-6);
258 Float_t intPtS = ptPara.Integral(fPtMin,fPtMax,1.e-6);
259#endif
6078e216 260 Float_t phiWgt=(fPhiMax-fPhiMin)/2./TMath::Pi(); //TR: should probably be done differently in case of anisotropic phi...
34f60c01 261 if (fAnalog == kAnalog) {
749070b6 262 fYWgt = intYS/fdNdy0;
263 fPtWgt = intPtS/intPt0;
264 fParentWeight = fYWgt*fPtWgt*phiWgt/fNpart;
265 } else {
266 fYWgt = intYS/fdNdy0;
267 fPtWgt = (fPtMax-fPtMin)/intPt0;
268 fParentWeight = fYWgt*fPtWgt*phiWgt/fNpart;
269 }
6078e216 270 //
271 // particle decay related initialization
00d6ce7d 272 fDecayer->SetForceDecay(fForceDecay);
fa480368 273 fDecayer->Init();
00d6ce7d 274
6078e216 275 //
fff02fee 276 AliGenMC::Init();
fe4da5cc 277}
278
279//____________________________________________________________
280void AliGenParam::Generate()
281{
d5615b85 282 //
283 // Generate 1 event (see Generate(Int_t ntimes) for details
6078e216 284 //
d5615b85 285 GenerateN(1);
286}
287//____________________________________________________________
288void AliGenParam::GenerateN(Int_t ntimes)
289{
290 //
291 // Generate ntimes*'npart' light and heavy mesons (J/Psi, upsilon or phi, Pion,
6078e216 292 // Kaons, Etas, Omegas) and Baryons (proton, antiprotons, neutrons and
293 // antineutrons in the the desired theta, phi and momentum windows;
294 // Gaussian smearing on the vertex is done if selected.
295 // The decay of heavy mesons is done using lujet,
296 // and the childern particle are tracked by GEANT
297 // However, light mesons are directly tracked by GEANT
298 // setting fForceDecay = nodecay (SetForceDecay(nodecay))
299 //
300 //
301 // Reinitialize decayer
2ad38d56 302 fDecayer->SetForceDecay(fForceDecay);
303 fDecayer->Init();
304
6078e216 305 //
28337bc1 306 Float_t polar[3]= {0,0,0}; // Polarisation of the parent particle (for GEANT tracking)
307 Float_t origin0[3]; // Origin of the generated parent particle (for GEANT tracking)
21391258 308 Float_t time0; // Time0 of the generated parent particle
28337bc1 309 Float_t pt, pl, ptot; // Transverse, logitudinal and total momenta of the parent particle
310 Float_t phi, theta; // Phi and theta spherical angles of the parent particle momentum
311 Float_t p[3], pc[3],
fff02fee 312 och[3]; // Momentum, polarisation and origin of the children particles from lujet
2d7a47be 313 Double_t ty, xmt;
fff02fee 314 Int_t nt, i, j;
fe4da5cc 315 Float_t wgtp, wgtch;
316 Double_t dummy;
09fd3ea2 317 static TClonesArray *particles;
fe4da5cc 318 //
fff02fee 319 if(!particles) particles = new TClonesArray("TParticle",1000);
8b31bfac 320
349be858 321 TDatabasePDG *pDataBase = TDatabasePDG::Instance();
fe4da5cc 322 //
323 Float_t random[6];
28337bc1 324
6078e216 325 // Calculating vertex position per event
fe4da5cc 326 for (j=0;j<3;j++) origin0[j]=fOrigin[j];
21391258 327 time0 = fTimeOrigin;
aee8290b 328 if(fVertexSmear==kPerEvent) {
d9ea0e3b 329 Vertex();
330 for (j=0;j<3;j++) origin0[j]=fVertex[j];
21391258 331 time0 = fTime;
fe4da5cc 332 }
d9ea0e3b 333
21aaa175 334 Int_t ipa=0;
fff02fee 335
6078e216 336 // Generating fNpart particles
7cf36cae 337 fNprimaries = 0;
338
d5615b85 339 Int_t nGen = fNpart*ntimes;
340 while (ipa<nGen) {
cc5d764c 341 while(1) {
6078e216 342 //
343 // particle type
65fb704d 344 Int_t iPart = fIpParaFunc(fRandom);
fa480368 345 fChildWeight=(fDecayer->GetPartialBranchingRatio(iPart))*fParentWeight;
675e9664 346 TParticlePDG *particle = pDataBase->GetParticle(iPart);
fa480368 347 Float_t am = particle->Mass();
8b31bfac 348
65fb704d 349 Rndm(random,2);
6078e216 350 //
351 // y
2d7a47be 352 ty = TMath::TanH(fYPara->GetRandom());
6078e216 353 //
354 // pT
34f60c01 355 if (fAnalog == kAnalog) {
fe4da5cc 356 pt=fPtPara->GetRandom();
357 wgtp=fParentWeight;
358 wgtch=fChildWeight;
359 } else {
360 pt=fPtMin+random[1]*(fPtMax-fPtMin);
361 Double_t ptd=pt;
362 wgtp=fParentWeight*fPtParaFunc(& ptd, &dummy);
363 wgtch=fChildWeight*fPtParaFunc(& ptd, &dummy);
364 }
365 xmt=sqrt(pt*pt+am*am);
2d7a47be 366 if (TMath::Abs(ty)==1.) {
367 ty=0.;
368 Fatal("AliGenParam",
369 "Division by 0: Please check you rapidity range !");
370 }
6078e216 371 //
372 // phi
373 // if(!ipa)
374 // phi=fEvPlane; //align first particle of each event with event plane
375 // else{
376 double v2 = fV2Para->Eval(pt);
377 fdNdPhi->SetParameter(0,v2);
378 fdNdPhi->SetParameter(1,fEvPlane);
379 phi=fdNdPhi->GetRandom();
380 // }
2d7a47be 381
60e55aee 382 pl=xmt*ty/sqrt((1.-ty)*(1.+ty));
fe4da5cc 383 theta=TMath::ATan2(pt,pl);
6078e216 384 // Cut on theta
fe4da5cc 385 if(theta<fThetaMin || theta>fThetaMax) continue;
386 ptot=TMath::Sqrt(pt*pt+pl*pl);
6078e216 387 // Cut on momentum
fe4da5cc 388 if(ptot<fPMin || ptot>fPMax) continue;
6078e216 389 //
fe4da5cc 390 p[0]=pt*TMath::Cos(phi);
391 p[1]=pt*TMath::Sin(phi);
392 p[2]=pl;
aee8290b 393 if(fVertexSmear==kPerTrack) {
65fb704d 394 Rndm(random,6);
fe4da5cc 395 for (j=0;j<3;j++) {
396 origin0[j]=
397 fOrigin[j]+fOsigma[j]*TMath::Cos(2*random[2*j]*TMath::Pi())*
398 TMath::Sqrt(-2*TMath::Log(random[2*j+1]));
399 }
21391258 400 Rndm(random,2);
401 time0 = fTimeOrigin + fOsigma[2]/TMath::Ccgs()*
402 TMath::Cos(2*random[0]*TMath::Pi())*
403 TMath::Sqrt(-2*TMath::Log(random[1]));
fe4da5cc 404 }
28337bc1 405
6078e216 406 // Looking at fForceDecay :
407 // if fForceDecay != none Primary particle decays using
408 // AliPythia and children are tracked by GEANT
409 //
410 // if fForceDecay == none Primary particle is tracked by GEANT
411 // (In the latest, make sure that GEANT actually does all the decays you want)
412 //
cd716030 413 Bool_t decayed = kFALSE;
414
fff02fee 415
34f60c01 416 if (fForceDecay != kNoDecay) {
6078e216 417 // Using lujet to decay particle
886b6f73 418 Float_t energy=TMath::Sqrt(ptot*ptot+am*am);
fa480368 419 TLorentzVector pmom(p[0], p[1], p[2], energy);
420 fDecayer->Decay(iPart,&pmom);
6078e216 421 //
422 // select decay particles
fa480368 423 Int_t np=fDecayer->ImportParticles(particles);
1242532d 424
425 // Selecting GeometryAcceptance for particles fPdgCodeParticleforAcceptanceCut;
426 if (fGeometryAcceptance)
427 if (!CheckAcceptanceGeometry(np,particles)) continue;
886b6f73 428 Int_t ncsel=0;
fff02fee 429 Int_t* pFlag = new Int_t[np];
430 Int_t* pParent = new Int_t[np];
431 Int_t* pSelected = new Int_t[np];
432 Int_t* trackIt = new Int_t[np];
433
434 for (i=0; i<np; i++) {
435 pFlag[i] = 0;
436 pSelected[i] = 0;
437 pParent[i] = -1;
438 }
439
6ba00c52 440 if (np >1) {
cd716030 441 decayed = kTRUE;
08bffa4d 442 TParticle* iparticle = 0;
fff02fee 443 Int_t ipF, ipL;
444 for (i = 1; i<np ; i++) {
445 trackIt[i] = 1;
6ba00c52 446 iparticle = (TParticle *) particles->At(i);
447 Int_t kf = iparticle->GetPdgCode();
fff02fee 448 Int_t ks = iparticle->GetStatusCode();
6078e216 449 // flagged particle
fff02fee 450
451 if (pFlag[i] == 1) {
fff02fee 452 ipF = iparticle->GetFirstDaughter();
453 ipL = iparticle->GetLastDaughter();
454 if (ipF > 0) for (j=ipF-1; j<ipL; j++) pFlag[j]=1;
455 continue;
456 }
457
6078e216 458 // flag decay products of particles with long life-time (c tau > .3 mum)
fff02fee 459
460 if (ks != 1) {
6078e216 461 // TParticlePDG *particle = pDataBase->GetParticle(kf);
fff02fee 462
463 Double_t lifeTime = fDecayer->GetLifetime(kf);
6078e216 464 // Double_t mass = particle->Mass();
465 // Double_t width = particle->Width();
47fc6bd5 466 if (lifeTime > (Double_t) fMaxLifeTime) {
fff02fee 467 ipF = iparticle->GetFirstDaughter();
468 ipL = iparticle->GetLastDaughter();
469 if (ipF > 0) for (j=ipF-1; j<ipL; j++) pFlag[j]=1;
470 } else{
471 trackIt[i] = 0;
472 pSelected[i] = 1;
473 }
474 } // ks==1 ?
6078e216 475 //
476 // children
47fc6bd5 477
18e09c20 478 if ((ChildSelected(TMath::Abs(kf)) || fForceDecay == kAll || fSelectAll) && trackIt[i])
6ba00c52 479 {
6ba00c52 480 if (fCutOnChild) {
fff02fee 481 pc[0]=iparticle->Px();
482 pc[1]=iparticle->Py();
483 pc[2]=iparticle->Pz();
484 Bool_t childok = KinematicSelection(iparticle, 1);
485 if(childok) {
486 pSelected[i] = 1;
6ba00c52 487 ncsel++;
488 } else {
489 ncsel=-1;
490 break;
491 } // child kine cuts
492 } else {
fff02fee 493 pSelected[i] = 1;
886b6f73 494 ncsel++;
6ba00c52 495 } // if child selection
496 } // select muon
497 } // decay particle loop
498 } // if decay products
499
886b6f73 500 Int_t iparent;
501 if ((fCutOnChild && ncsel >0) || !fCutOnChild){
502 ipa++;
6078e216 503 //
504 // Parent
cd716030 505
506
21391258 507 PushTrack(0, -1, iPart, p, origin0, polar, time0, kPPrimary, nt, wgtp, ((decayed)? 11 : 1));
fff02fee 508 pParent[0] = nt;
a99cf51f 509 KeepTrack(nt);
7cf36cae 510 fNprimaries++;
511
6078e216 512 //
513 // Decay Products
514 //
fff02fee 515 for (i = 1; i < np; i++) {
516 if (pSelected[i]) {
517 TParticle* iparticle = (TParticle *) particles->At(i);
ea79897e 518 Int_t kf = iparticle->GetPdgCode();
cd716030 519 Int_t ksc = iparticle->GetStatusCode();
ea79897e 520 Int_t jpa = iparticle->GetFirstMother()-1;
fff02fee 521
774ceaaf 522 och[0] = origin0[0]+iparticle->Vx();
523 och[1] = origin0[1]+iparticle->Vy();
524 och[2] = origin0[2]+iparticle->Vz();
fff02fee 525 pc[0] = iparticle->Px();
526 pc[1] = iparticle->Py();
527 pc[2] = iparticle->Pz();
528
ea79897e 529 if (jpa > -1) {
530 iparent = pParent[jpa];
fff02fee 531 } else {
532 iparent = -1;
533 }
1242532d 534
ea79897e 535 PushTrack(fTrackIt * trackIt[i], iparent, kf,
fff02fee 536 pc, och, polar,
21391258 537 time0 + iparticle->T(), kPDecay, nt, wgtch, ksc);
fff02fee 538 pParent[i] = nt;
a99cf51f 539 KeepTrack(nt);
7cf36cae 540 fNprimaries++;
8b31bfac 541 } // Selected
542 } // Particle loop
28337bc1 543 } // Decays by Lujet
8b31bfac 544 particles->Clear();
fff02fee 545 if (pFlag) delete[] pFlag;
546 if (pParent) delete[] pParent;
547 if (pSelected) delete[] pSelected;
8b31bfac 548 if (trackIt) delete[] trackIt;
21aaa175 549 } // kinematic selection
28337bc1 550 else // nodecay option, so parent will be tracked by GEANT (pions, kaons, eta, omegas, baryons)
551 {
5d12ce38 552 gAlice->GetMCApp()->
21391258 553 PushTrack(fTrackIt,-1,iPart,p,origin0,polar,time0,kPPrimary,nt,wgtp, 1);
28337bc1 554 ipa++;
7cf36cae 555 fNprimaries++;
28337bc1 556 }
fe4da5cc 557 break;
21aaa175 558 } // while
fe4da5cc 559 } // event loop
7cf36cae 560
a99cf51f 561 SetHighWaterMark(nt);
7cf36cae 562
563 AliGenEventHeader* header = new AliGenEventHeader("PARAM");
564 header->SetPrimaryVertex(fVertex);
21391258 565 header->SetInteractionTime(fTime);
7cf36cae 566 header->SetNProduced(fNprimaries);
567 AddHeader(header);
fe4da5cc 568}
2ad38d56 569//____________________________________________________________________________________
570Float_t AliGenParam::GetRelativeArea(Float_t ptMin, Float_t ptMax, Float_t yMin, Float_t yMax, Float_t phiMin, Float_t phiMax)
571{
6078e216 572 //
573 // Normalisation for selected kinematic region
574 //
faf00b0d 575#if ROOT_VERSION_CODE < ROOT_VERSION(5,99,0)
2ad38d56 576 Float_t ratio =
389a9124 577 fPtPara->Integral(ptMin,ptMax,(Double_t *)0,1.e-6) / fPtPara->Integral( fPtPara->GetXmin(), fPtPara->GetXmax(),(Double_t *)0,1.e-6) *
578 fYPara->Integral(yMin,yMax,(Double_t *)0,1.e-6)/fYPara->Integral(fYPara->GetXmin(),fYPara->GetXmax(),(Double_t *)0,1.e-6) *
2ad38d56 579 (phiMax-phiMin)/360.;
faf00b0d 580#else
581 Float_t ratio =
582 fPtPara->Integral(ptMin,ptMax,1.e-6) / fPtPara->Integral( fPtPara->GetXmin(), fPtPara->GetXmax(),1.e-6) *
583 fYPara->Integral(yMin,yMax,1.e-6)/fYPara->Integral(fYPara->GetXmin(),fYPara->GetXmax(),1.e-6) *
584 (phiMax-phiMin)/360.;
585#endif
2ad38d56 586 return TMath::Abs(ratio);
587}
588
589//____________________________________________________________________________________
fe4da5cc 590
dc1d768c 591void AliGenParam::Draw( const char * /*opt*/)
5bd39445 592{
593 //
594 // Draw the pT and y Distributions
595 //
596 TCanvas *c0 = new TCanvas("c0","Canvas 0",400,10,600,700);
597 c0->Divide(2,1);
598 c0->cd(1);
599 fPtPara->Draw();
600 fPtPara->GetHistogram()->SetXTitle("p_{T} (GeV)");
601 c0->cd(2);
602 fYPara->Draw();
603 fYPara->GetHistogram()->SetXTitle("y");
604}
605
f87cfe57 606
fe4da5cc 607
608