// //
////////////////////////////////////////////////////////////////////////////
+#include <iostream>
+#include <string>
+
#include "TUHKMgen.h"
#ifndef DATABASE_PDG
#include "DatabasePDG.h"
#include <TParticle.h>
#include <TClonesArray.h>
#include <TMCProcess.h>
-#include "TDatabasePDG.h"
-#include "TSystem.h"
+#include <TDatabasePDG.h>
+#include <TSystem.h>
#include "AliGenUHKM.h"
#include "AliRun.h"
#include "AliGenHijingEventHeader.h"
#include "AliLog.h"
-#include <iostream>
-#include <string>
using namespace std;
fHydjetParams.fIenglu=0;
fHydjetParams.fIanglu=0;
*/
- strcpy(fParticleFilename, Form("%s/TUHKMgen/UHKM/particles.data", gSystem->Getenv("ALICE_ROOT")));
- strcpy(fDecayFilename, Form("%s/TUHKMgen/UHKM/tabledecay.txt", gSystem->Getenv("ALICE_ROOT")));
+ strncpy(fParticleFilename, Form("%s/TUHKMgen/UHKM/particles.data", gSystem->Getenv("ALICE_ROOT")), 255);
+ strncpy(fDecayFilename, Form("%s/TUHKMgen/UHKM/tabledecay.txt", gSystem->Getenv("ALICE_ROOT")), 255);
for(Int_t i=0; i<500; i++) {
fStableFlagPDG[i] = 0;
fStableFlagStatus[i] = kFALSE;
fHydjetParams.fIanglu=0;
*/
- strcpy(fParticleFilename, Form("%s/TUHKMgen/UHKM/particles.data", gSystem->Getenv("ALICE_ROOT")));
- strcpy(fDecayFilename, Form("%s/TUHKMgen/UHKM/tabledecay.txt", gSystem->Getenv("ALICE_ROOT")));
+ strncpy(fParticleFilename, Form("%s/TUHKMgen/UHKM/particles.data", gSystem->Getenv("ALICE_ROOT")), 255);
+ strncpy(fDecayFilename, Form("%s/TUHKMgen/UHKM/tabledecay.txt", gSystem->Getenv("ALICE_ROOT")), 255);
for(Int_t i=0; i<500; i++) {
fStableFlagPDG[i] = 0;
fStableFlagStatus[i] = kFALSE;
Float_t polar[3] = {0,0,0};
Float_t origin[3] = {0,0,0};
Float_t origin0[3] = {0,0,0};
+ Float_t time0 = 0.;
+ Float_t p[3];
Float_t v[3];
- Float_t mass, energy;
+ Float_t mass=0.0, energy=0.0;
Vertex();
for(Int_t j=0; j<3; j++) origin0[j] = fVertex[j];
+ time0 = fTime;
- fTrials = 0;
-
- Int_t nt = 0;
-
+ // Generate the event and import particles
fUHKMgen->GenerateEvent();
- fTrials++;
-
fUHKMgen->ImportParticles(&fParticles,"All");
-
Int_t np = fParticles.GetEntriesFast();
+ Int_t nt = 0;
+ // Handle the IDs of particles on the stack
Int_t* idsOnStack = new Int_t[np];
Int_t* newPos = new Int_t[np];
- for(Int_t i=0; i<np; i++) newPos[i] = i;
+ for(Int_t i=0; i<np; i++) {
+ newPos[i] = i;
+ idsOnStack[i] = -1;
+ }
+
+
+ // Generate a random phi used to rotate the whole event
+ Double_t eventRotation = gRandom->Rndm()*TMath::Pi();
- //_________ Loop for particle selection
+ TParticle *iparticle;
+ Double_t partMomPhi=0.0, partPt=0.0;
+ Double_t partVtxPhi=0.0, partVtxR=0.0;
+ //_________ Loop for particles in the stack
for(Int_t i=0; i<np; i++) {
- TParticle *iparticle = (TParticle*)fParticles.At(i);
+ iparticle = (TParticle*)fParticles.At(i);
Int_t kf = iparticle->GetPdgCode();
Bool_t hasMother = (iparticle->GetFirstMother() >= 0);
Bool_t hasDaughter = (iparticle->GetNDaughters() > 0);
// It will not be tracked
// Add it only once with coordinates not
// smeared with primary vertex position
- Float_t p[3] = {p[0] = iparticle->Px(),
- p[1] = iparticle->Py(),
- p[2] = iparticle->Pz()};
+
+ // rotate the direction of the particle
+ partMomPhi = TMath::ATan2(iparticle->Py(), iparticle->Px());
+ partPt = TMath::Hypot(iparticle->Px(), iparticle->Py());
+ p[0] = partPt*TMath::Cos(partMomPhi+eventRotation);
+ p[1] = partPt*TMath::Sin(partMomPhi+eventRotation);
+ p[2] = iparticle->Pz();
+
mass = TDatabasePDG::Instance()->GetParticle(kf)->Mass();
energy = sqrt(mass*mass + p[0]*p[0] + p[1]*p[1] + p[2]*p[2]);
- v[0] = iparticle->Vx();
- v[1] = iparticle->Vy();
+
+ // rotate the freezeout point
+ partVtxPhi = TMath::ATan2(iparticle->Vy(), iparticle->Vx());
+ partVtxR = TMath::Hypot(iparticle->Vx(), iparticle->Vy());
+ v[0] = partVtxR*TMath::Cos(partVtxPhi + eventRotation);
+ v[1] = partVtxR*TMath::Cos(partVtxPhi + eventRotation);
v[2] = iparticle->Vz();
- Float_t time = iparticle->T();
+ Float_t time = iparticle->T();
+
Int_t imo = -1;
if(hasMother) {
- imo = iparticle->GetFirstMother(); //index of mother particle in fParticles
+ imo = iparticle->GetFirstMother(); //index of mother particle in fParticles
} // if has mother
Bool_t trackFlag = kFALSE; // tFlag is kFALSE --> do not track the particle
PushTrack(trackFlag, (imo>=0 ? idsOnStack[imo] : imo), kf,
- p[0], p[1], p[2], energy,
- v[0], v[1], v[2], time,
- polar[0], polar[1], polar[2],
- (hasMother ? kPDecay : kPNoProcess), nt);
+ p[0], p[1], p[2], energy,
+ v[0], v[1], v[2], time,
+ polar[0], polar[1], polar[2],
+ (hasMother ? kPDecay : kPNoProcess), nt);
idsOnStack[i] = nt;
fNprimaries++;
// this one will not be tracked
// Second time with event-wide c0ordinates and vertex smearing
// this one will be tracked
-
- Float_t p[3] = {p[0] = iparticle->Px(),
- p[1] = iparticle->Py(),
- p[2] = iparticle->Pz()};
- mass = TDatabasePDG::Instance()->GetParticle(kf)->Mass();
+
+ // rotate the direction of the particle
+ partMomPhi = TMath::ATan2(iparticle->Py(), iparticle->Px());
+ partPt = TMath::Hypot(iparticle->Px(), iparticle->Py());
+ p[0] = partPt*TMath::Cos(partMomPhi+eventRotation);
+ p[1] = partPt*TMath::Sin(partMomPhi+eventRotation);
+ p[2] = iparticle->Pz();
+
energy = sqrt(mass*mass + p[0]*p[0] + p[1]*p[1] + p[2]*p[2]);
- v[0] = iparticle->Vx();
- v[1] = iparticle->Vy();
- v[2] = iparticle->Vz();
+ // rotate the freezeout point
+ partVtxPhi = TMath::ATan2(iparticle->Vy(), iparticle->Vx());
+ partVtxR = TMath::Hypot(iparticle->Vx(), iparticle->Vy());
+ v[0] = partVtxR*TMath::Cos(partVtxPhi + eventRotation);
+ v[1] = partVtxR*TMath::Cos(partVtxPhi + eventRotation);
+ v[2] = iparticle->Vz();
+
Int_t type = iparticle->GetStatusCode(); // 1-from jet / 0-from hydro
Int_t coeffT=1;
if(type==1) coeffT=-1; //to separate particles from jets
-
+
Int_t imo = -1;
if(hasMother) {
- imo = iparticle->GetFirstMother();
+ imo = iparticle->GetFirstMother();
} // if has mother
Bool_t trackFlag = kFALSE; // tFlag = kFALSE --> do not track this one, its for femtoscopy
PushTrack(trackFlag, (imo>=0 ? idsOnStack[imo] : imo), kf,
- p[0], p[1], p[2], energy,
- v[0], v[1], v[2], (iparticle->T())*coeffT,
- polar[0], polar[1], polar[2],
- hasMother ? kPDecay:kPNoProcess, nt);
+ p[0], p[1], p[2], energy,
+ v[0], v[1], v[2], (iparticle->T())*coeffT, // freeze-out time is negative if the particle comes from jet
+ polar[0], polar[1], polar[2],
+ hasMother ? kPDecay:kPNoProcess, nt);
idsOnStack[i] = nt;
fNprimaries++;
KeepTrack(nt);
-
+
origin[0] = origin0[0]+v[0];
origin[1] = origin0[1]+v[1];
origin[2] = origin0[2]+v[2];
+ Float_t time = time0+iparticle->T();
imo = nt;
- trackFlag = kTRUE; // tFlag = kTRUE --> track this one
-
+ trackFlag = fTrackIt; // Track this particle, unless otherwise specified by fTrackIt
+
PushTrack(trackFlag, imo, kf,
- p[0], p[1], p[2], energy,
- origin[0], origin[1], origin[2], iparticle->T(),
- polar[0], polar[1], polar[2],
- hasMother ? kPDecay:kPNoProcess, nt);
+ p[0], p[1], p[2], energy,
+ origin[0], origin[1], origin[2], time,
+ polar[0], polar[1], polar[2],
+ hasMother ? kPDecay:kPNoProcess, nt);
fNprimaries++;
KeepTrack(nt);
}
}
-
+
SetHighWaterMark(fNprimaries);
TArrayF eventVertex;
eventVertex[0] = origin0[0];
eventVertex[1] = origin0[1];
eventVertex[2] = origin0[2];
+ Float_t eventTime = time0;
// Builds the event header, to be called after each event
AliGenEventHeader* header = new AliGenHijingEventHeader("UHKM");
+ Double_t b = 0.;
+ Double_t npart = 0;
+ Double_t nbin = 0;
+ fUHKMgen->GetCentrality(b, npart, nbin);
+ printf("********** %13.3f %13.3f %13.3f \n", b, npart, nbin);
+
+
((AliGenHijingEventHeader*) header)->SetNProduced(fNprimaries);
((AliGenHijingEventHeader*) header)->SetPrimaryVertex(eventVertex);
- ((AliGenHijingEventHeader*) header)->SetImpactParameter(0.0);
+ ((AliGenHijingEventHeader*) header)->SetInteractionTime(eventTime);
+ ((AliGenHijingEventHeader*) header)->SetImpactParameter(b);
((AliGenHijingEventHeader*) header)->SetTotalEnergy(0.0);
((AliGenHijingEventHeader*) header)->SetHardScatters(0);
- ((AliGenHijingEventHeader*) header)->SetParticipants(0, 0);
- ((AliGenHijingEventHeader*) header)->SetCollisions(0, 0, 0, 0);
+ ((AliGenHijingEventHeader*) header)->SetParticipants(Int_t(npart), 0);
+ ((AliGenHijingEventHeader*) header)->SetCollisions(Int_t(nbin), 0, 0, 0);
((AliGenHijingEventHeader*) header)->SetSpectators(0, 0, 0, 0);
((AliGenHijingEventHeader*) header)->SetReactionPlaneAngle(0);//evrot);
header->SetPrimaryVertex(fVertex);
+ header->SetInteractionTime(fTime);
AddHeader(header);
fCollisionGeometry = (AliGenHijingEventHeader*) header;
- delete idsOnStack;
-
+ delete [] idsOnStack;
+ delete [] newPos;
}
void AliGenUHKM::Copy(TObject &) const
<< TDatabasePDG::Instance()->GetParticle(uhkmTestParticle->GetPDG())->Width() << endl;
}
} // end for
-}
\ No newline at end of file
+}
git://git.uio.no / u / mrichter / AliRoot.git / blobdiff