/*
$Log$
+Revision 1.32 2001/01/26 19:55:51 hristov
+Major upgrade of AliRoot code
+
+Revision 1.31 2001/01/17 10:54:31 hristov
+Better protection against FPE
+
+Revision 1.30 2000/12/18 08:55:35 morsch
+Make AliPythia dependent generartors work with new scheme of random number generation
+
+Revision 1.29 2000/12/04 11:22:03 morsch
+Init of sRandom as in 1.15
+
+Revision 1.28 2000/12/02 11:41:39 morsch
+Use SetRandom() to initialize random number generator in constructor.
+
+Revision 1.27 2000/11/30 20:29:02 morsch
+Initialise static variable sRandom in constructor: sRandom = fRandom;
+
+Revision 1.26 2000/11/30 07:12:50 alibrary
+Introducing new Rndm and QA classes
+
+Revision 1.25 2000/10/18 19:11:27 hristov
+Division by zero fixed
+
+Revision 1.24 2000/09/18 10:41:35 morsch
+Add possibility to use nuclear structure functions from PDF library V8.
+
+Revision 1.23 2000/09/14 14:05:40 morsch
+dito
+
+Revision 1.22 2000/09/14 14:02:22 morsch
+- Correct conversion from mm to cm when passing particle vertex to MC.
+- Correct handling of fForceDecay == all.
+
+Revision 1.21 2000/09/12 14:14:55 morsch
+Call fDecayer->ForceDecay() at the beginning of Generate().
+
+Revision 1.20 2000/09/06 14:29:33 morsch
+Use AliPythia for event generation an AliDecayPythia for decays.
+Correct handling of "nodecay" option
+
+Revision 1.19 2000/07/11 18:24:56 fca
+Coding convention corrections + few minor bug fixes
+
+Revision 1.18 2000/06/30 12:40:34 morsch
+Pythia takes care of vertex smearing. Correct conversion from Pythia units (mm) to
+Geant units (cm).
+
+Revision 1.17 2000/06/09 20:34:07 morsch
+All coding rule violations except RS3 corrected
+
+Revision 1.16 2000/05/15 15:04:20 morsch
+The full event is written for fNtrack = -1
+Coding rule violations corrected.
+
Revision 1.15 2000/04/26 10:14:24 morsch
Particles array has one entry more than pythia particle list. Upper bound of
particle loop changed to np-1 (R. Guernane, AM)
Introduction of the Copyright and cvs Log
*/
-#include "AliGenerator.h"
#include "AliGenPythia.h"
+#include "AliDecayerPythia.h"
#include "AliRun.h"
#include "AliPythia.h"
-#include <TDirectory.h>
-#include <TFile.h>
-#include <TTree.h>
-#include <stdlib.h>
-#include <AliPythia.h>
+#include "AliPDG.h"
#include <TParticle.h>
-//#include <GParticle.h>
+#include <TSystem.h>
+
ClassImp(AliGenPythia)
AliGenPythia::AliGenPythia()
:AliGenerator()
{
-// Constructor
+// Default Constructor
+ fDecayer = new AliDecayerPythia();
}
AliGenPythia::AliGenPythia(Int_t npart)
// structure function GRVHO
//
fXsection = 0.;
+ fNucA1=0;
+ fNucA2=0;
fParentSelect.Set(5);
fChildSelect.Set(5);
for (Int_t i=0; i<5; i++) fParentSelect[i]=fChildSelect[i]=0;
SetForceDecay();
SetPtHard();
SetEnergyCMS();
+ fDecayer = new AliDecayerPythia();
+ // Set random number generator
+ sRandom=fRandom;
+}
+
+AliGenPythia::AliGenPythia(const AliGenPythia & Pythia)
+{
+// copy constructor
}
AliGenPythia::~AliGenPythia()
void AliGenPythia::Init()
{
// Initialisation
- SetMC(new AliPythia());
+ SetMC(AliPythia::Instance());
fPythia=(AliPythia*) fgMCEvGen;
//
fParentWeight=1./Float_t(fNpart);
//
// Forward Paramters to the AliPythia object
- fPythia->DefineParticles();
+ // gSystem->Exec("ln -s $ALICE_ROOT/data/Decay.table fort.1");
+ // fPythia->Pyupda(2,1);
+ // gSystem->Exec("rm fort.1");
+
+ fDecayer->SetForceDecay(fForceDecay);
+ fDecayer->Init();
+
+
fPythia->SetCKIN(3,fPtHardMin);
fPythia->SetCKIN(4,fPtHardMax);
+ if (fNucA1 > 0 && fNucA2 > 0) fPythia->SetNuclei(fNucA1, fNucA2);
fPythia->ProcInit(fProcess,fEnergyCMS,fStrucFunc);
- fPythia->ForceDecay(fForceDecay);
- fPythia->Lulist(0);
- fPythia->Pystat(2);
+
+ // fPythia->Pylist(0);
+ // fPythia->Pystat(2);
// Parent and Children Selection
switch (fProcess)
{
case dielectron:
case b_jpsi_dielectron:
case b_psip_dielectron:
- fChildSelect[0]=11;
+ fChildSelect[0]=kElectron;
break;
case semimuonic:
case dimuon:
case b_psip_dimuon:
case pitomu:
case katomu:
- fChildSelect[0]=13;
+ fChildSelect[0]=kMuonMinus;
break;
+ case hadronicD:
+ fChildSelect[0]=kPiPlus;
+ fChildSelect[1]=kKPlus;
+ break;
case all:
case nodecay:
- break;
+ break;
}
}
void AliGenPythia::Generate()
{
// Generate one event
+ fDecayer->ForceDecay();
Float_t polar[3] = {0,0,0};
Float_t origin[3]= {0,0,0};
Float_t originP[3]= {0,0,0};
Float_t origin0[3]= {0,0,0};
- Float_t p[3], pP[4], random[6];
+ Float_t p[3], pP[4];
+// Float_t random[6];
static TClonesArray *particles;
// converts from mm/c to s
const Float_t kconv=0.001/2.999792458e8;
fTrials=0;
for (j=0;j<3;j++) origin0[j]=fOrigin[j];
- if(fVertexSmear==perEvent) {
- gMC->Rndm(random,6);
+ if(fVertexSmear==kPerEvent) {
+ fPythia->SetMSTP(151,1);
for (j=0;j<3;j++) {
- origin0[j]+=fOsigma[j]*TMath::Cos(2*random[2*j]*TMath::Pi())*
- TMath::Sqrt(-2*TMath::Log(random[2*j+1]));
- fPythia->SetMSTP(151,0);
+ fPythia->SetPARP(151+j, fOsigma[j]/10.);
}
- } else if (fVertexSmear==perTrack) {
+ } else if (fVertexSmear==kPerTrack) {
fPythia->SetMSTP(151,0);
- for (j=0;j<3;j++) {
- fPythia->SetPARP(151+j, fOsigma[j]*10.);
- }
}
+
while(1)
{
fPythia->Pyevnt();
- fPythia->Lulist(1);
+// fPythia->Pylist(1);
fTrials++;
fPythia->ImportParticles(particles,"All");
Int_t np = particles->GetEntriesFast();
kf = CheckPDGCode(iparticle->GetPdgCode());
if (ks==21) continue;
- fChildWeight=(fPythia->GetBraPart(kf))*fParentWeight;
+ fChildWeight=(fDecayer->GetPartialBranchingRatio(kf))*fParentWeight;
//
// Parent
if (ParentSelected(TMath::Abs(kf))) {
pP[3] = massP;
gAlice->SetTrack(0,-1,-1,
pP,originP,polar,
- 0,"Hard Scat.",ntP,fParentWeight);
+ 0,kPPrimary,ntP,fParentWeight);
+// 0,"Hard Scat.",ntP,fParentWeight);
gAlice->KeepTrack(ntP);
}
nc++;
p[0]=iparticle->Px();
p[1]=iparticle->Py();
p[2]=iparticle->Pz();
- origin[0]=origin0[0]+iparticle->Vx()/10;
- origin[1]=origin0[1]+iparticle->Vy()/10;
- origin[2]=origin0[2]+iparticle->Vz()/10;
+ origin[0]=origin0[0]+iparticle->Vx()/10.;
+ origin[1]=origin0[1]+iparticle->Vy()/10.;
+ origin[2]=origin0[2]+iparticle->Vz()/10.;
Int_t ifch=iparticle->GetFirstDaughter();
Int_t ilch=iparticle->GetLastDaughter();
- if (ifch !=0 && ilch !=0) {
- gAlice->SetTrack(0,ntP,kf,
+
+ if ((ifch !=0 && ilch !=0) || fForceDecay == nodecay) {
+ Int_t trackit=0;
+ if (fForceDecay == nodecay) trackit = 1;
+ gAlice->SetTrack(trackit,ntP,kf,
p,origin,polar,
- 0,"Primary",nt,fParentWeight);
+ 0,kPPrimary,nt,fParentWeight);
gAlice->KeepTrack(nt);
Int_t iparent = nt;
//
// Children
-
- for (j=ifch; j<=ilch; j++)
- {
- TParticle * ichild =
+ if (fForceDecay != nodecay) {
+ for (j=ifch; j<=ilch; j++)
+ {
+ TParticle * ichild =
(TParticle *) particles->At(j-1);
- kf = CheckPDGCode(ichild->GetPdgCode());
+ kf = CheckPDGCode(ichild->GetPdgCode());
//
//
- if (ChildSelected(TMath::Abs(kf))) {
- origin[0]=ichild->Vx();
- origin[1]=ichild->Vy();
- origin[2]=ichild->Vz();
+ if (ChildSelected(TMath::Abs(kf))) {
+ origin[0]=origin0[0]+ichild->Vx()/10.;
+ origin[1]=origin0[1]+ichild->Vy()/10.;
+ origin[2]=origin0[2]+ichild->Vz()/10.;
p[0]=ichild->Px();
p[1]=ichild->Py();
p[2]=ichild->Pz();
Float_t tof=kconv*ichild->T();
gAlice->SetTrack(fTrackIt, iparent, kf,
p,origin,polar,
- tof,"Decay",nt,fChildWeight);
+ tof,kPDecay,nt,fChildWeight);
gAlice->KeepTrack(nt);
- } // select child
- } // child loop
+ } // select child
+ } // child loop
+ }
}
} // kinematic selection
} // select particle
p[0]=iparticle->Px();
p[1]=iparticle->Py();
p[2]=iparticle->Pz();
- origin[0]=origin0[0]+iparticle->Vx()/10;
- origin[1]=origin0[1]+iparticle->Vy()/10;
- origin[2]=origin0[2]+iparticle->Vz()/10;
+ origin[0]=origin0[0]+iparticle->Vx()/10.;
+ origin[1]=origin0[1]+iparticle->Vy()/10.;
+ origin[2]=origin0[2]+iparticle->Vz()/10.;
Float_t tof=kconv*iparticle->T();
gAlice->SetTrack(fTrackIt,-1,kf,p,origin,polar,
- tof,"Primary",nt);
+ tof,kPPrimary,nt);
gAlice->KeepTrack(nt);
} // select particle
} // particle loop
if (jev >= fNpart || fNpart == -1) {
fKineBias=Float_t(fNpart)/Float_t(fTrials);
printf("\n Trials: %i %i %i\n",fTrials, fNpart, jev);
+// Print x-section summary
+ fPythia->Pystat(1);
break;
}
}
Bool_t AliGenPythia::ChildSelected(Int_t ip)
{
// True if particle is in list of decay products to be selected
+ if (fForceDecay == all) return kTRUE;
+
for (Int_t i=0; i<5; i++)
{
if (fChildSelect[i]==ip) return kTRUE;
//
// rapidity cut
- Float_t y = 0.5*TMath::Log((e+pz)/(e-pz));
- if (y > fYMax || y < fYMin)
- {
+ if ( (e-pz)<=0 || (e+pz)<=0 ) {
+ return kFALSE;
+ }
+ else {
+ Float_t y = 0.5*TMath::Log((e+pz)/(e-pz));
+ if (y > fYMax || y < fYMin)
+ {
// printf("\n failed y cut %f %f %f \n",y,fYMin,fYMax);
- return kFALSE;
+ return kFALSE;
+ }
}
//
// phi cut
- Float_t phi=Float_t(TMath::ATan2(Double_t(py),Double_t(px)))+TMath::Pi();
+ Float_t phi=Float_t(TMath::ATan2(Double_t(py),Double_t(px)));
if (phi > fPhiMax || phi < fPhiMin)
{
// printf("\n failed phi cut %f %f %f \n",phi,fPhiMin,fPhiMax);
{
// Adjust the weights after generation of all events
//
- TClonesArray *partArray = gAlice->Particles();
TParticle *part;
Int_t ntrack=gAlice->GetNtrack();
for (Int_t i=0; i<ntrack; i++) {
- part= (TParticle*) partArray->UncheckedAt(i);
+ part= gAlice->Particle(i);
part->SetWeight(part->GetWeight()*fKineBias);
}
}
//
// If the particle is in a diffractive state, then take action accordingly
switch (pdgcode) {
+ case 91:
+ return 92;
case 110:
//rho_diff0 -- difficult to translate, return rho0
return 113;
//non diffractive state -- return code unchanged
return pdgcode;
}
-
-
-
-
+
+void AliGenPythia::SetNuclei(Int_t a1, Int_t a2)
+{
+// Treat protons as inside nuclei with mass numbers a1 and a2
+ fNucA1 = a1;
+ fNucA2 = a2;
+}
+
+
+AliGenPythia& AliGenPythia::operator=(const AliGenPythia& rhs)
+{
+// Assignment operator
+ return *this;
+}
+#ifdef never
+void AliGenPythia::Streamer(TBuffer &R__b)
+{
+ // Stream an object of class AliGenPythia.
+
+ if (R__b.IsReading()) {
+ Version_t R__v = R__b.ReadVersion(); if (R__v) { }
+ AliGenerator::Streamer(R__b);
+ R__b >> (Int_t&)fProcess;
+ R__b >> (Int_t&)fStrucFunc;
+ R__b >> (Int_t&)fForceDecay;
+ R__b >> fEnergyCMS;
+ R__b >> fKineBias;
+ R__b >> fTrials;
+ fParentSelect.Streamer(R__b);
+ fChildSelect.Streamer(R__b);
+ R__b >> fXsection;
+// (AliPythia::Instance())->Streamer(R__b);
+ R__b >> fPtHardMin;
+ R__b >> fPtHardMax;
+// if (fDecayer) fDecayer->Streamer(R__b);
+ } else {
+ R__b.WriteVersion(AliGenPythia::IsA());
+ AliGenerator::Streamer(R__b);
+ R__b << (Int_t)fProcess;
+ R__b << (Int_t)fStrucFunc;
+ R__b << (Int_t)fForceDecay;
+ R__b << fEnergyCMS;
+ R__b << fKineBias;
+ R__b << fTrials;
+ fParentSelect.Streamer(R__b);
+ fChildSelect.Streamer(R__b);
+ R__b << fXsection;
+// R__b << fPythia;
+ R__b << fPtHardMin;
+ R__b << fPtHardMax;
+ // fDecayer->Streamer(R__b);
+ }
+}
+#endif
git://git.uio.no / u / mrichter / AliRoot.git / blobdiff