/*
$Log$
+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
*/
#include "AliGenPythia.h"
+#include "AliDecayerPythia.h"
#include "AliRun.h"
#include "AliPythia.h"
+#include "AliPDG.h"
#include <TParticle.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)
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};
fTrials=0;
for (j=0;j<3;j++) origin0[j]=fOrigin[j];
- if(fVertexSmear==perEvent) {
+ if(fVertexSmear==kPerEvent) {
fPythia->SetMSTP(151,1);
for (j=0;j<3;j++) {
fPythia->SetPARP(151+j, fOsigma[j]/10.);
}
- } else if (fVertexSmear==perTrack) {
+ } else if (fVertexSmear==kPerTrack) {
fPythia->SetMSTP(151,0);
}
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]=origin0[0]+ichild->Vx()*10.;
- origin[1]=origin0[1]+ichild->Vy()*10.;
- origin[2]=origin0[2]+ichild->Vz()*10.;
+ 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) {
+ 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;
+ }
}
//
//
// 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
-
-
+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);
+ }
+}