#include "AliPythia.h"
#include "AliPythiaRndm.h"
#include "AliRun.h"
+#include "AliStack.h"
+#include "AliRunLoader.h"
+#include "AliMC.h"
- ClassImp(AliGenPythia)
+ClassImp(AliGenPythia)
AliGenPythia::AliGenPythia()
:AliGenMC()
SetGammaPhiRange();
SetGammaEtaRange();
SetPtKick();
+ SetQuench();
+
fSetNuclei = kFALSE;
if (!AliPythiaRndm::GetPythiaRandom())
AliPythiaRndm::SetPythiaRandom(GetRandom());
SetGammaPhiRange();
SetGammaEtaRange();
SetJetReconstructionMode();
+ SetQuench();
SetPtKick();
// Options determining what to keep in the stack (Heavy flavour generation)
fStackFillOpt = kFlavorSelection; // Keep particle with selected flavor
}
AliGenPythia::AliGenPythia(const AliGenPythia & Pythia)
+ :AliGenMC(Pythia)
{
// copy constructor
Pythia.Copy(*this);
SetMC(AliPythia::Instance());
fPythia=(AliPythia*) fMCEvGen;
+ printf("Pythia Init %p", fPythia);
+
//
fParentWeight=1./Float_t(fNpart);
//
Float_t polar[3] = {0,0,0};
Float_t origin[3] = {0,0,0};
- Float_t p[3];
+ Float_t p[4];
// converts from mm/c to s
const Float_t kconv=0.001/2.999792458e8;
//
// event loop
while(1)
{
+ if (fQuench) {
+ fPythia->SetMSTJ(1, 0);
+ }
+ printf("Pythia pyevnt %p", fPythia);
fPythia->Pyevnt();
- if (gAlice->GetEvNumber()>=fDebugEventFirst &&
- gAlice->GetEvNumber()<=fDebugEventLast) fPythia->Pylist(1);
- fTrials++;
+
+ if (fQuench) {
+ fPythia->Quench();
+ fPythia->SetMSTJ(1, 1);
+ fPythia->Pyexec();
+ }
+
+ if (gAlice) {
+ if (gAlice->GetEvNumber()>=fDebugEventFirst &&
+ gAlice->GetEvNumber()<=fDebugEventLast) fPythia->Pylist(1);
+ }
+ fTrials++;
+ printf("Pythia import t %p %p ", fPythia, fParticles);
fPythia->ImportParticles(fParticles,"All");
Boost();
p[0] = iparticle->Px();
p[1] = iparticle->Py();
p[2] = iparticle->Pz();
+ p[3] = iparticle->Energy();
+
origin[0] = fOrigin[0]+iparticle->Vx()/10.;
origin[1] = fOrigin[1]+iparticle->Vy()/10.;
origin[2] = fOrigin[2]+iparticle->Vz()/10.;
Float_t tof = kconv*iparticle->T();
Int_t ipa = iparticle->GetFirstMother()-1;
Int_t iparent = (ipa > -1) ? pParent[ipa] : -1;
- PushTrack(fTrackIt*trackIt[i] ,
- iparent, kf, p, origin, polar, tof, kPPrimary, nt, 1., ks);
+
+ PushTrack(fTrackIt*trackIt[i], iparent, kf,
+ p[0], p[1], p[2], p[3],
+ origin[0], origin[1], origin[2], tof,
+ polar[0], polar[1], polar[2],
+ kPPrimary, nt, 1., ks);
pParent[i] = nt;
KeepTrack(nt);
} // PushTrack loop
p[0] = iparticle->Px();
p[1] = iparticle->Py();
p[2] = iparticle->Pz();
+ p[3] = iparticle->Energy();
+
origin[0] = fOrigin[0]+iparticle->Vx()/10.;
origin[1] = fOrigin[1]+iparticle->Vy()/10.;
origin[2] = fOrigin[2]+iparticle->Vz()/10.;
Float_t tof=kconv*iparticle->T();
- PushTrack(fTrackIt*trackIt, iparent, kf, p, origin, polar,
- tof, kPPrimary, nt, 1., ks);
+
+ PushTrack(fTrackIt*trackIt, iparent, kf,
+ p[0], p[1], p[2], p[3],
+ origin[0], origin[1], origin[2], tof,
+ polar[0], polar[1], polar[2],
+ kPPrimary, nt, 1., ks);
KeepTrack(nt);
pParent[i] = nt;
} // select particle
{
// Adjust the weights after generation of all events
//
- TParticle *part;
- Int_t ntrack=gAlice->GetNtrack();
- for (Int_t i=0; i<ntrack; i++) {
- part= gAlice->Particle(i);
- part->SetWeight(part->GetWeight()*fKineBias);
+ if (gAlice) {
+ TParticle *part;
+ Int_t ntrack=gAlice->GetMCApp()->GetNtrack();
+ for (Int_t i=0; i<ntrack; i++) {
+ part= gAlice->GetMCApp()->Particle(i);
+ part->SetWeight(part->GetWeight()*fKineBias);
+ }
}
}
jets[3][i]);
}
}
- gAlice->SetGenEventHeader(fHeader);
+ if (gAlice) gAlice->SetGenEventHeader(fHeader);
}
AliGenPythia& AliGenPythia::operator=(const AliGenPythia& rhs)
{
// Assignment operator
+ rhs.Copy(*this);
return *this;
}
// Load event into Pythia Common Block
//
-
- Int_t npart = (Int_t) (gAlice->TreeK())->GetEntries();
+ AliRunLoader* rl = AliRunLoader::GetRunLoader();
+ Int_t npart = (rl->Stack())-> GetNprimary();
(fPythia->GetPyjets())->N = npart;
for (Int_t part = 0; part < npart; part++) {
- TParticle *MPart = gAlice->Particle(part);
+ TParticle *MPart = (rl->Stack())->Particle(part);
Int_t kf = MPart->GetPdgCode();
Int_t ks = MPart->GetStatusCode();
Float_t px = MPart->Px();
Float_t py = MPart->Py();
Float_t pz = MPart->Pz();
Float_t e = MPart->Energy();
- Float_t p = TMath::Sqrt(px * px + py * py + pz * pz);
- Float_t m = TMath::Sqrt(e * e - p * p);
+ Float_t m = MPart->GetCalcMass();
(fPythia->GetPyjets())->P[0][part] = px;
}
}
-void AliGenPythia::RecJetsUA1(Float_t eCellMin, Float_t eCellSeed, Float_t eMin, Float_t rMax,
- Int_t& njets, Float_t jets [4][50])
+void AliGenPythia::RecJetsUA1(Int_t& njets, Float_t jets [4][50])
{
//
// Calls the Pythia jet finding algorithm to find jets in the current event
Float_t pz = (fPythia->GetPyjets())->P[2][n+i];
Float_t e = (fPythia->GetPyjets())->P[3][n+i];
Float_t pt = TMath::Sqrt(px * px + py * py);
- Float_t phi = TMath::ATan2(py,px);
+ Float_t phi = TMath::Pi() + TMath::ATan2(-py, -px);
Float_t theta = TMath::ATan2(pt,pz);
Float_t et = e * TMath::Sin(theta);
Float_t eta = -TMath::Log(TMath::Tan(theta / 2.));